module CNPhenologyMod

#include "shr_assert.h"

  !-----------------------------------------------------------------------
  ! !MODULE: CNPhenologyMod
  !
  ! !DESCRIPTION:
  ! Module holding routines used in phenology model for coupled carbon
  ! nitrogen code.
  !
  ! !USES:
  use shr_kind_mod                    , only : r8 => shr_kind_r8
  use shr_log_mod                     , only : errMsg => shr_log_errMsg
  use shr_sys_mod                     , only : shr_sys_flush
  use decompMod                       , only : bounds_type
  use clm_varpar                      , only : maxveg, nlevdecomp_full
  use clm_varctl                      , only : iulog, use_cndv
  use clm_varcon                      , only : tfrz
  use abortutils                      , only : endrun
  use CanopyStateType                 , only : canopystate_type
  use CNDVType                        , only : dgvs_type
  use CNVegstateType                  , only : cnveg_state_type
  use CNVegCarbonStateType            , only : cnveg_carbonstate_type
  use CNVegCarbonFluxType             , only : cnveg_carbonflux_type
  use CNVegnitrogenstateType          , only : cnveg_nitrogenstate_type
  use CNVegnitrogenfluxType           , only : cnveg_nitrogenflux_type
  use CropType                        , only : crop_type
  use pftconMod                       , only : pftcon
  use SoilStateType                   , only : soilstate_type
  use TemperatureType                 , only : temperature_type
  use WaterDiagnosticBulkType         , only : waterdiagnosticbulk_type
  use Wateratm2lndBulkType            , only : wateratm2lndbulk_type
  use initVerticalMod                 , only : find_soil_layer_containing_depth
  use ColumnType                      , only : col
  use GridcellType                    , only : grc                
  use PatchType                       , only : patch   
  use atm2lndType                     , only : atm2lnd_type             
  !
  implicit none
  private
  !
  ! !PUBLIC MEMBER FUNCTIONS:
  public :: readParams           ! Read parameters
  public :: CNPhenologyreadNML   ! Read namelist
  public :: CNPhenologyInit      ! Initialization
  public :: CNPhenology          ! Update
  !
  ! !PRIVATE DATA MEMBERS:
  type, private :: params_type
     real(r8) :: crit_dayl       ! critical day length for senescence
     real(r8) :: ndays_on     	 ! number of days to complete leaf onset
     real(r8) :: ndays_off	 ! number of days to complete leaf offset
     real(r8) :: fstor2tran      ! fraction of storage to move to transfer for each onset
     real(r8) :: crit_onset_fdd  ! critical number of freezing days to set gdd counter
     real(r8) :: crit_onset_swi  ! critical number of days > soilpsi_on for onset
     real(r8) :: soilpsi_on      ! critical soil water potential for leaf onset
     real(r8) :: crit_offset_fdd ! critical number of freezing days to initiate offset
     real(r8) :: crit_offset_swi ! critical number of water stress days to initiate offset
     real(r8) :: soilpsi_off     ! critical soil water potential for leaf offset
     real(r8) :: lwtop   	 ! live wood turnover proportion (annual fraction)
     real(r8) :: phenology_soil_depth ! soil depth used for measuring states for phenology triggers
  end type params_type

  type(params_type) :: params_inst

  real(r8) :: dt                            ! radiation time step delta t (seconds)
  real(r8) :: fracday                       ! dtime as a fraction of day
  real(r8) :: crit_dayl                     ! critical daylength for offset (seconds)
  real(r8) :: ndays_on                      ! number of days to complete onset
  real(r8) :: ndays_off                     ! number of days to complete offset
  real(r8) :: fstor2tran                    ! fraction of storage to move to transfer on each onset
  real(r8) :: crit_onset_fdd                ! critical number of freezing days
  real(r8) :: crit_onset_swi                ! water stress days for offset trigger
  real(r8) :: soilpsi_on                    ! water potential for onset trigger (MPa)
  real(r8) :: crit_offset_fdd               ! critical number of freezing degree days to trigger offset
  real(r8) :: crit_offset_swi               ! water stress days for offset trigger
  real(r8) :: soilpsi_off                   ! water potential for offset trigger (MPa)
  real(r8) :: lwtop                         ! live wood turnover proportion (annual fraction)
  integer  :: phenology_soil_layer          ! soil layer used for measuring states for phenology triggers

  ! CropPhenology variables and constants
  real(r8) :: p1d, p1v                      ! photoperiod factor constants for crop vernalization
  real(r8) :: hti                           ! cold hardening index threshold for vernalization
  real(r8) :: tbase                         ! base temperature for vernalization

  integer, parameter :: NOT_Planted   = 999 ! If not planted   yet in year
  integer, parameter :: NOT_Harvested = 999 ! If not harvested yet in year
  integer, parameter :: inNH       = 1      ! Northern Hemisphere
  integer, parameter :: inSH       = 2      ! Southern Hemisphere
  integer, pointer   :: inhemi(:)           ! Hemisphere that patch is in 

  integer, allocatable :: minplantjday(:,:) ! minimum planting julian day
  integer, allocatable :: maxplantjday(:,:) ! maximum planting julian day
  integer              :: jdayyrstart(inSH) ! julian day of start of year

  real(r8), private :: initial_seed_at_planting = 3._r8 ! Initial seed at planting

  character(len=*), parameter, private :: sourcefile = &
       __FILE__
  !-----------------------------------------------------------------------

contains

  !-----------------------------------------------------------------------
  subroutine CNPhenologyReadNML( NLFilename )
    !
    ! !DESCRIPTION:
    ! Read the namelist for CNPhenology
    !
    ! !USES:
    use fileutils      , only : getavu, relavu, opnfil
    use shr_nl_mod     , only : shr_nl_find_group_name
    use spmdMod        , only : masterproc, mpicom
    use shr_mpi_mod    , only : shr_mpi_bcast
    use clm_varctl     , only : iulog
    !
    ! !ARGUMENTS:
    character(len=*), intent(in) :: NLFilename ! Namelist filename
    !
    ! !LOCAL VARIABLES:
    integer :: ierr                 ! error code
    integer :: unitn                ! unit for namelist file

    character(len=*), parameter :: subname = 'CNPhenologyReadNML'
    character(len=*), parameter :: nmlname = 'cnphenology'
    !-----------------------------------------------------------------------
    namelist /cnphenology/ initial_seed_at_planting

    ! Initialize options to default values, in case they are not specified in
    ! the namelist

    if (masterproc) then
       unitn = getavu()
       write(iulog,*) 'Read in '//nmlname//'  namelist'
       call opnfil (NLFilename, unitn, 'F')
       call shr_nl_find_group_name(unitn, nmlname, status=ierr)
       if (ierr == 0) then
          read(unitn, nml=cnphenology, iostat=ierr)
          if (ierr /= 0) then
             call endrun(msg="ERROR reading "//nmlname//"namelist"//errmsg(sourcefile, __LINE__))
          end if
       else
          call endrun(msg="ERROR could NOT find "//nmlname//"namelist"//errmsg(sourcefile, __LINE__))
       end if
       call relavu( unitn )
    end if

    call shr_mpi_bcast (initial_seed_at_planting, mpicom)

    if (masterproc) then
       write(iulog,*) ' '
       write(iulog,*) nmlname//' settings:'
       write(iulog,nml=cnphenology)
       write(iulog,*) ' '
    end if


    !-----------------------------------------------------------------------
    
  end subroutine CNPhenologyReadNML

  !-----------------------------------------------------------------------
  subroutine readParams ( ncid )
    !
    ! !DESCRIPTION:
    !
    ! !USES:
    use ncdio_pio    , only: file_desc_t
    use paramUtilMod , only : readNcdioScalar

    ! !ARGUMENTS:
    implicit none
    type(file_desc_t),intent(inout) :: ncid   ! pio netCDF file id
    !
    ! !LOCAL VARIABLES:
    character(len=*), parameter  :: subname = 'readParams_CNPhenology'
    !-----------------------------------------------------------------------

    call readNcdioScalar(ncid, 'crit_dayl', subname, params_inst%crit_dayl)
    call readNcdioScalar(ncid, 'ndays_on', subname, params_inst%ndays_on)
    call readNcdioScalar(ncid, 'ndays_off', subname, params_inst%ndays_off)
    call readNcdioScalar(ncid, 'fstor2tran', subname, params_inst%fstor2tran)
    call readNcdioScalar(ncid, 'crit_onset_fdd', subname, params_inst%crit_onset_fdd)
    call readNcdioScalar(ncid, 'crit_onset_swi', subname, params_inst%crit_onset_swi)
    call readNcdioScalar(ncid, 'soilpsi_on', subname, params_inst%soilpsi_on)
    call readNcdioScalar(ncid, 'crit_offset_fdd', subname, params_inst%crit_offset_fdd)
    call readNcdioScalar(ncid, 'crit_offset_swi', subname, params_inst%crit_offset_swi)
    call readNcdioScalar(ncid, 'soilpsi_off', subname, params_inst%soilpsi_off)
    call readNcdioScalar(ncid, 'lwtop_ann', subname, params_inst%lwtop)
    call readNcdioScalar(ncid, 'phenology_soil_depth', subname, params_inst%phenology_soil_depth)

  end subroutine readParams

  !-----------------------------------------------------------------------
  subroutine CNPhenology (bounds, num_soilc, filter_soilc, num_soilp, &
       filter_soilp, num_pcropp, filter_pcropp, &
       doalb, waterdiagnosticbulk_inst, wateratm2lndbulk_inst, temperature_inst, atm2lnd_inst, crop_inst, &
       canopystate_inst, soilstate_inst, dgvs_inst, &
       cnveg_state_inst, cnveg_carbonstate_inst, cnveg_carbonflux_inst,    &
       cnveg_nitrogenstate_inst, cnveg_nitrogenflux_inst, &
       c13_cnveg_carbonstate_inst, c14_cnveg_carbonstate_inst, &
       leaf_prof_patch, froot_prof_patch, phase)
    ! !USES:
    use CNSharedParamsMod, only: use_fun
    !
    ! !DESCRIPTION:
    ! Dynamic phenology routine for coupled carbon-nitrogen code (CN)
    ! 1. grass phenology
    !
    ! !ARGUMENTS:
    type(bounds_type)              , intent(in)    :: bounds
    integer                        , intent(in)    :: num_soilc       ! number of soil columns in filter
    integer                        , intent(in)    :: filter_soilc(:) ! filter for soil columns
    integer                        , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                        , intent(in)    :: filter_soilp(:) ! filter for soil patches
    integer                        , intent(in)    :: num_pcropp      ! number of prog. crop patches in filter
    integer                        , intent(in)    :: filter_pcropp(:)! filter for prognostic crop patches
    logical                        , intent(in)    :: doalb           ! true if time for sfc albedo calc
    type(waterdiagnosticbulk_type)          , intent(in)    :: waterdiagnosticbulk_inst
    type(wateratm2lndbulk_type)          , intent(in)    :: wateratm2lndbulk_inst
    type(temperature_type)         , intent(inout) :: temperature_inst
    type(atm2lnd_type)             , intent(in)    :: atm2lnd_inst
    type(crop_type)                , intent(inout) :: crop_inst
    type(canopystate_type)         , intent(in)    :: canopystate_inst
    type(soilstate_type)           , intent(in)    :: soilstate_inst
    type(dgvs_type)                , intent(inout) :: dgvs_inst
    type(cnveg_state_type)         , intent(inout) :: cnveg_state_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: cnveg_carbonstate_inst
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenstate_type) , intent(inout) :: cnveg_nitrogenstate_inst
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: c13_cnveg_carbonstate_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: c14_cnveg_carbonstate_inst
    real(r8)                       , intent(in)    :: leaf_prof_patch(bounds%begp:,1:)
    real(r8)                       , intent(in)    :: froot_prof_patch(bounds%begp:,1:)
    integer                        , intent(in)    :: phase
    !-----------------------------------------------------------------------

    SHR_ASSERT_ALL_FL((ubound(leaf_prof_patch)   == (/bounds%endp,nlevdecomp_full/)), sourcefile, __LINE__)
    SHR_ASSERT_ALL_FL((ubound(froot_prof_patch)  == (/bounds%endp,nlevdecomp_full/)), sourcefile, __LINE__)

    ! each of the following phenology type routines includes a filter
    ! to operate only on the relevant patches


    if ( phase == 1 ) then
       call CNPhenologyClimate(num_soilp, filter_soilp, num_pcropp, filter_pcropp, &
            temperature_inst, cnveg_state_inst, crop_inst)
   
       call CNEvergreenPhenology(num_soilp, filter_soilp, &
            cnveg_state_inst, cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       call CNSeasonDecidPhenology(num_soilp, filter_soilp, &
            temperature_inst, cnveg_state_inst, dgvs_inst, &
            cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       call CNStressDecidPhenology(num_soilp, filter_soilp,   &
            soilstate_inst, temperature_inst, atm2lnd_inst, wateratm2lndbulk_inst, cnveg_state_inst, &
            cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       if (doalb .and. num_pcropp > 0 ) then
          call CropPhenology(num_pcropp, filter_pcropp, &
               waterdiagnosticbulk_inst, temperature_inst, crop_inst, canopystate_inst, cnveg_state_inst, &
               cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst, &
               c13_cnveg_carbonstate_inst, c14_cnveg_carbonstate_inst)
       end if
    else if ( phase == 2 ) then
       ! the same onset and offset routines are called regardless of
       ! phenology type - they depend only on onset_flag, offset_flag, bglfr, and bgtr

       call CNOnsetGrowth(num_soilp, filter_soilp, &
            cnveg_state_inst, &
            cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       call CNOffsetLitterfall(num_soilp, filter_soilp, &
            cnveg_state_inst, cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       call CNBackgroundLitterfall(num_soilp, filter_soilp, &
            cnveg_state_inst, cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
   
       call CNLivewoodTurnover(num_soilp, filter_soilp, &
            cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       call CNCropHarvestToProductPools(bounds, num_soilp, filter_soilp, num_soilc, filter_soilc, &
            cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)

       ! gather all patch-level litterfall fluxes to the column for litter C and N inputs

       call CNLitterToColumn(bounds, num_soilc, filter_soilc, &
            cnveg_state_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst, &
            leaf_prof_patch(bounds%begp:bounds%endp,1:nlevdecomp_full), & 
            froot_prof_patch(bounds%begp:bounds%endp,1:nlevdecomp_full))
    else
       call endrun( 'bad phase' )
    end if

  end subroutine CNPhenology

  !-----------------------------------------------------------------------
  subroutine CNPhenologyInit(bounds)
    !
    ! !DESCRIPTION:
    ! Initialization of CNPhenology. Must be called after time-manager is
    ! initialized, and after pftcon file is read in.
    !
    ! !USES:
    use clm_time_manager, only: get_step_size_real
    use clm_varctl      , only: use_crop
    use clm_varcon      , only: secspday
    !
    ! !ARGUMENTS:
    type(bounds_type), intent(in) :: bounds  
    !------------------------------------------------------------------------

    !
    ! Get time-step and what fraction of a day it is
    !
    dt      = get_step_size_real()
    fracday = dt/secspday

    ! set constants for CNSeasonDecidPhenology 
    ! (critical daylength from Biome-BGC, v4.1.2)
    crit_dayl=params_inst%crit_dayl

    ! Set constants for CNSeasonDecidPhenology and CNStressDecidPhenology
    ndays_on=params_inst%ndays_on
    ndays_off=params_inst%ndays_off

    ! set transfer parameters
    fstor2tran=params_inst%fstor2tran

    call find_soil_layer_containing_depth( &
         depth = params_inst%phenology_soil_depth, &
         layer = phenology_soil_layer)

    ! -----------------------------------------
    ! Constants for CNStressDecidPhenology
    ! -----------------------------------------

    ! onset parameters
    crit_onset_fdd=params_inst%crit_onset_fdd
    ! critical onset gdd now being calculated as a function of annual
    ! average 2m temp.
    ! crit_onset_gdd = 150.0 ! c3 grass value
    ! crit_onset_gdd = 1000.0   ! c4 grass value
    crit_onset_swi=params_inst%crit_onset_swi
    soilpsi_on=params_inst%soilpsi_on

    ! offset parameters
    crit_offset_fdd=params_inst%crit_offset_fdd
    crit_offset_swi=params_inst%crit_offset_swi
    soilpsi_off=params_inst%soilpsi_off    

    ! -----------------------------------------
    ! Constants for CNLivewoodTurnover
    ! -----------------------------------------

    ! set the global parameter for livewood turnover rate
    ! define as an annual fraction (0.7), and convert to fraction per second
    lwtop=params_inst%lwtop/31536000.0_r8 !annual fraction converted to per second

    ! -----------------------------------------
    ! Call any subroutine specific initialization routines
    ! -----------------------------------------

    if ( use_crop ) call CropPhenologyInit(bounds)

  end subroutine CNPhenologyInit

  !-----------------------------------------------------------------------
  subroutine CNPhenologyClimate (num_soilp, filter_soilp, num_pcropp, filter_pcropp, &
       temperature_inst, cnveg_state_inst, crop_inst)
    !
    ! !DESCRIPTION:
    ! For coupled carbon-nitrogen code (CN).
    !
    ! !USES:
    use clm_time_manager , only : get_days_per_year
    use clm_time_manager , only : get_curr_date, is_first_step
    !
    ! !ARGUMENTS:
    integer                , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                , intent(in)    :: filter_soilp(:) ! filter for soil patches
    integer                , intent(in)    :: num_pcropp      ! number of prognostic crops in filter
    integer                , intent(in)    :: filter_pcropp(:)! filter for prognostic crop patches
    type(temperature_type) , intent(inout) :: temperature_inst
    type(cnveg_state_type) , intent(inout) :: cnveg_state_inst
    type(crop_type)        , intent(inout) :: crop_inst
    !
    ! !LOCAL VARIABLES:
    integer  :: p       ! indices
    integer  :: fp      ! lake filter patch index
    real(r8) :: dayspyr ! days per year (days)
    integer  :: kyr     ! current year
    integer  :: kmo     ! month of year  (1, ..., 12)
    integer  :: kda     ! day of month   (1, ..., 31)
    integer  :: mcsec   ! seconds of day (0, ..., seconds/day)
    real(r8), parameter :: yravg   = 20.0_r8      ! length of years to average for gdd
    real(r8), parameter :: yravgm1 = yravg-1.0_r8 ! minus 1 of above
    !-----------------------------------------------------------------------

    associate(                                                & 
         nyrs_crop_active => crop_inst%nyrs_crop_active_patch,   & ! InOut:  [integer (:)  ]  number of years this crop patch has been active
         
         t_ref2m        => temperature_inst%t_ref2m_patch ,   & ! Input:  [real(r8) (:) ]  2m air temperature (K)
         gdd0           => temperature_inst%gdd0_patch    ,   & ! Output: [real(r8) (:) ]  growing deg. days base 0 deg C (ddays)            
         gdd8           => temperature_inst%gdd8_patch    ,   & ! Output: [real(r8) (:) ]     "     "    "    "   8  "  "    "               
         gdd10          => temperature_inst%gdd10_patch   ,   & ! Output: [real(r8) (:) ]     "     "    "    "  10  "  "    "               
         gdd020         => temperature_inst%gdd020_patch  ,   & ! Output: [real(r8) (:) ]  20-yr mean of gdd0 (ddays)                        
         gdd820         => temperature_inst%gdd820_patch  ,   & ! Output: [real(r8) (:) ]  20-yr mean of gdd8 (ddays)                        
         gdd1020        => temperature_inst%gdd1020_patch ,   & ! Output: [real(r8) (:) ]  20-yr mean of gdd10 (ddays)                       
         
         tempavg_t2m    => cnveg_state_inst%tempavg_t2m_patch & ! Output: [real(r8) (:) ]  temp. avg 2m air temperature (K)                  
         )

      ! set time steps
      
      dayspyr = get_days_per_year()

      do fp = 1,num_soilp
         p = filter_soilp(fp)
         tempavg_t2m(p) = tempavg_t2m(p) + t_ref2m(p) * (fracday/dayspyr)
      end do

      !
      ! The following crop related steps are done here rather than CropPhenology
      ! so that they will be completed each time-step rather than with doalb.
      !
      ! The following lines come from ibis's climate.f + stats.f
      ! gdd SUMMATIONS ARE RELATIVE TO THE PLANTING DATE (see subr. updateAccFlds)

      if (num_pcropp > 0) then
         ! get time-related info
         call get_curr_date(kyr, kmo, kda, mcsec)
      end if

      do fp = 1,num_pcropp
         p = filter_pcropp(fp)
         if (kmo == 1 .and. kda == 1 .and. nyrs_crop_active(p) == 0) then ! YR 1:
            gdd020(p)  = 0._r8                             ! set gdd..20 variables to 0
            gdd820(p)  = 0._r8                             ! and crops will not be planted
            gdd1020(p) = 0._r8
         end if
         if (kmo == 1 .and. kda == 1 .and. mcsec == 0) then        ! <-- END of EVERY YR:
            if (nyrs_crop_active(p) == 1) then                     ! <-- END of YR 1
               gdd020(p)  = gdd0(p)                                ! <-- END of YR 1
               gdd820(p)  = gdd8(p)                                ! <-- END of YR 1
               gdd1020(p) = gdd10(p)                               ! <-- END of YR 1
            end if                                                 ! <-- END of YR 1
            gdd020(p)  = (yravgm1* gdd020(p)  + gdd0(p))  / yravg  ! gdd..20 must be long term avgs
            gdd820(p)  = (yravgm1* gdd820(p)  + gdd8(p))  / yravg  ! so ignore results for yrs 1 & 2
            gdd1020(p) = (yravgm1* gdd1020(p) + gdd10(p)) / yravg 
         end if
      end do

    end associate

  end subroutine CNPhenologyClimate

  !-----------------------------------------------------------------------
  subroutine CNEvergreenPhenology (num_soilp, filter_soilp , &
       cnveg_state_inst, cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)   
       ! cnveg_state_inst) 
    !
    ! !DESCRIPTION:
    ! For coupled carbon-nitrogen code (CN).
    !
    ! !USES:
    use clm_varcon       , only : secspday
    use clm_time_manager , only : get_days_per_year
    use clm_varctl       , only : CN_evergreen_phenology_opt   
    !
    ! !ARGUMENTS:
    integer           , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer           , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(cnveg_state_type), intent(inout) :: cnveg_state_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: cnveg_carbonstate_inst    
    type(cnveg_nitrogenstate_type) , intent(inout) :: cnveg_nitrogenstate_inst  
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst     
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst   
    !
    ! !LOCAL VARIABLES:
    real(r8):: dayspyr                    ! Days per year
    integer :: p                          ! indices
    integer :: fp                         ! lake filter patch index
    
    real(r8):: tranr 				      
    real(r8):: t1                         ! temporary variable 
    !-----------------------------------------------------------------------

    associate(                                        & 
         ivt        => patch%itype                    , & ! Input:  [integer  (:) ]  patch vegetation type                                

         evergreen  => pftcon%evergreen             , & ! Input:  binary flag for evergreen leaf habit (0 or 1)     
         leaf_long  => pftcon%leaf_long             , & ! Input:  leaf longevity (yrs)  
         
         woody                               =>    pftcon%woody                                                         , & ! Input:  binary flag for woody lifeform (1=woody, 0=not woody)     
         
         leafc_storage                       =>    cnveg_carbonstate_inst%leafc_storage_patch                           , & ! Input:  [real(r8) (:)]  (gC/m2) leaf C storage                             
   		 frootc_storage                      =>    cnveg_carbonstate_inst%frootc_storage_patch                          , & ! Input:  [real(r8) (:)]  (gC/m2) fine root C storage                         
         livestemc_storage                   =>    cnveg_carbonstate_inst%livestemc_storage_patch                       , & ! Input:  [real(r8) (:)]  (gC/m2) live stem C storage                         
         deadstemc_storage                   =>    cnveg_carbonstate_inst%deadstemc_storage_patch                       , & ! Input:  [real(r8) (:)]  (gC/m2) dead stem C storage                         
   		 livecrootc_storage                  =>    cnveg_carbonstate_inst%livecrootc_storage_patch                      , & ! Input:  [real(r8) (:)]  (gC/m2) live coarse root C storage                  
   		 deadcrootc_storage                  =>    cnveg_carbonstate_inst%deadcrootc_storage_patch                      , & ! Input:  [real(r8) (:)]  (gC/m2) dead coarse root C storage                  
   		 gresp_storage                       =>    cnveg_carbonstate_inst%gresp_storage_patch                           , & ! Input:  [real(r8) (:)]  (gC/m2) growth respiration storage   
   		 leafc_xfer                          =>    cnveg_carbonstate_inst%leafc_xfer_patch                              , & ! InOut:  [real(r8) (:)]  (gC/m2) leaf C transfer                            
   		 frootc_xfer                         =>    cnveg_carbonstate_inst%frootc_xfer_patch                             , & ! InOut:  [real(r8) (:)]  (gC/m2) fine root C transfer                       
   		 livestemc_xfer                      =>    cnveg_carbonstate_inst%livestemc_xfer_patch                          , & ! InOut:  [real(r8) (:)]  (gC/m2) live stem C transfer                       
   		 deadstemc_xfer                      =>    cnveg_carbonstate_inst%deadstemc_xfer_patch                          , & ! InOut:  [real(r8) (:)]  (gC/m2) dead stem C transfer                       
   		 livecrootc_xfer                     =>    cnveg_carbonstate_inst%livecrootc_xfer_patch                         , & ! InOut:  [real(r8) (:)]  (gC/m2) live coarse root C transfer                
   		 deadcrootc_xfer                     =>    cnveg_carbonstate_inst%deadcrootc_xfer_patch                         , & ! InOut:  [real(r8) (:)]  (gC/m2) dead coarse root C transfer   
   		                
   		 leafn_storage                       =>    cnveg_nitrogenstate_inst%leafn_storage_patch                         , & ! Input:  [real(r8) (:)]  (gN/m2) leaf N storage                              
   		 frootn_storage                      =>    cnveg_nitrogenstate_inst%frootn_storage_patch                        , & ! Input:  [real(r8) (:)]  (gN/m2) fine root N storage                         
   		 livestemn_storage                   =>    cnveg_nitrogenstate_inst%livestemn_storage_patch                     , & ! Input:  [real(r8) (:)]  (gN/m2) live stem N storage                         
   		 deadstemn_storage                   =>    cnveg_nitrogenstate_inst%deadstemn_storage_patch                     , & ! Input:  [real(r8) (:)]  (gN/m2) dead stem N storage                         
   		 livecrootn_storage                  =>    cnveg_nitrogenstate_inst%livecrootn_storage_patch                    , & ! Input:  [real(r8) (:)]  (gN/m2) live coarse root N storage                  
   		 deadcrootn_storage                  =>    cnveg_nitrogenstate_inst%deadcrootn_storage_patch                    , & ! Input:  [real(r8) (:)]  (gN/m2) dead coarse root N storage               
   		 leafn_xfer                          =>    cnveg_nitrogenstate_inst%leafn_xfer_patch                            , & ! InOut:  [real(r8) (:)]  (gN/m2) leaf N transfer                            
   		 frootn_xfer                         =>    cnveg_nitrogenstate_inst%frootn_xfer_patch                           , & ! InOut:  [real(r8) (:)]  (gN/m2) fine root N transfer                       
   		 livestemn_xfer                      =>    cnveg_nitrogenstate_inst%livestemn_xfer_patch                        , & ! InOut:  [real(r8) (:)]  (gN/m2) live stem N transfer                       
   		 deadstemn_xfer                      =>    cnveg_nitrogenstate_inst%deadstemn_xfer_patch                        , & ! InOut:  [real(r8) (:)]  (gN/m2) dead stem N transfer                       
   		 livecrootn_xfer                     =>    cnveg_nitrogenstate_inst%livecrootn_xfer_patch                       , & ! InOut:  [real(r8) (:)]  (gN/m2) live coarse root N transfer                
   		 deadcrootn_xfer                     =>    cnveg_nitrogenstate_inst%deadcrootn_xfer_patch                       , & ! InOut:  [real(r8) (:)]  (gN/m2) dead coarse root N transfer     
   		                 
   		 leafc_storage_to_xfer               =>    cnveg_carbonflux_inst%leafc_storage_to_xfer_patch                    , & ! InOut:  [real(r8) (:)]                                                     
   		 frootc_storage_to_xfer              =>    cnveg_carbonflux_inst%frootc_storage_to_xfer_patch                   , & ! InOut:  [real(r8) (:)]                                                     
   		 livestemc_storage_to_xfer           =>    cnveg_carbonflux_inst%livestemc_storage_to_xfer_patch                , & ! InOut:  [real(r8) (:)]                                                     
   		 deadstemc_storage_to_xfer           =>    cnveg_carbonflux_inst%deadstemc_storage_to_xfer_patch                , & ! InOut:  [real(r8) (:)]                                                     
   		 livecrootc_storage_to_xfer          =>    cnveg_carbonflux_inst%livecrootc_storage_to_xfer_patch               , & ! InOut:  [real(r8) (:)]                                                     
   		 deadcrootc_storage_to_xfer          =>    cnveg_carbonflux_inst%deadcrootc_storage_to_xfer_patch               , & ! InOut:  [real(r8) (:)]                                                     
   		 gresp_storage_to_xfer               =>    cnveg_carbonflux_inst%gresp_storage_to_xfer_patch                    , & ! InOut:  [real(r8) (:)]  
   		 leafc_xfer_to_leafc                 =>    cnveg_carbonflux_inst%leafc_xfer_to_leafc_patch                      , & ! InOut:  [real(r8) (:)]                                                    
   		 frootc_xfer_to_frootc               =>    cnveg_carbonflux_inst%frootc_xfer_to_frootc_patch                    , & ! InOut:  [real(r8) (:)]                                                    
   		 livestemc_xfer_to_livestemc         =>    cnveg_carbonflux_inst%livestemc_xfer_to_livestemc_patch              , & ! InOut:  [real(r8) (:)]                                                    
   		 deadstemc_xfer_to_deadstemc         =>    cnveg_carbonflux_inst%deadstemc_xfer_to_deadstemc_patch              , & ! InOut:  [real(r8) (:)]                                                    
   		 livecrootc_xfer_to_livecrootc       =>    cnveg_carbonflux_inst%livecrootc_xfer_to_livecrootc_patch            , & ! InOut:  [real(r8) (:)]                                                    
   		 deadcrootc_xfer_to_deadcrootc       =>    cnveg_carbonflux_inst%deadcrootc_xfer_to_deadcrootc_patch            , & ! InOut:  [real(r8) (:)]   
   		                                                    
   		 leafn_storage_to_xfer               =>    cnveg_nitrogenflux_inst%leafn_storage_to_xfer_patch                  , & ! InOut:  [real(r8) (:)]                                                     
   		 frootn_storage_to_xfer              =>    cnveg_nitrogenflux_inst%frootn_storage_to_xfer_patch                 , & ! InOut:  [real(r8) (:)]                                                     
   		 livestemn_storage_to_xfer           =>    cnveg_nitrogenflux_inst%livestemn_storage_to_xfer_patch              , & ! InOut:  [real(r8) (:)]                                                     
   		 deadstemn_storage_to_xfer           =>    cnveg_nitrogenflux_inst%deadstemn_storage_to_xfer_patch              , & ! InOut:  [real(r8) (:)]                                                     
   		 livecrootn_storage_to_xfer          =>    cnveg_nitrogenflux_inst%livecrootn_storage_to_xfer_patch             , & ! InOut:  [real(r8) (:)]   
   		 deadcrootn_storage_to_xfer          =>    cnveg_nitrogenflux_inst%deadcrootn_storage_to_xfer_patch             , & ! InOut:  [real(r8) (:)]                                                     
   		 leafn_xfer_to_leafn                 =>    cnveg_nitrogenflux_inst%leafn_xfer_to_leafn_patch                    , & ! InOut:  [real(r8) (:)]                                                    
   		 frootn_xfer_to_frootn               =>    cnveg_nitrogenflux_inst%frootn_xfer_to_frootn_patch                  , & ! InOut:  [real(r8) (:)]                                                    
   		 livestemn_xfer_to_livestemn         =>    cnveg_nitrogenflux_inst%livestemn_xfer_to_livestemn_patch            , & ! InOut:  [real(r8) (:)]                                                    
   		 deadstemn_xfer_to_deadstemn         =>    cnveg_nitrogenflux_inst%deadstemn_xfer_to_deadstemn_patch            , & ! InOut:  [real(r8) (:)]                                                    
   		 livecrootn_xfer_to_livecrootn       =>    cnveg_nitrogenflux_inst%livecrootn_xfer_to_livecrootn_patch          , & ! InOut:  [real(r8) (:)]                                                    
   		 deadcrootn_xfer_to_deadcrootn       =>    cnveg_nitrogenflux_inst%deadcrootn_xfer_to_deadcrootn_patch          , & ! InOut:  [real(r8) (:)]     
   		           
         bglfr      => cnveg_state_inst%bglfr_patch , & ! Output: [real(r8) (:) ]  background litterfall rate (1/s)                  
         bgtr       => cnveg_state_inst%bgtr_patch  , & ! Output: [real(r8) (:) ]  background transfer growth rate (1/s)             
         lgsf       => cnveg_state_inst%lgsf_patch    & ! Output: [real(r8) (:) ]  long growing season factor [0-1]                  
         )

      dayspyr   = get_days_per_year()

      do fp = 1,num_soilp
         p = filter_soilp(fp)
         if (evergreen(ivt(p)) == 1._r8) then
            bglfr(p) = 1._r8/(leaf_long(ivt(p)) * dayspyr * secspday)
            bgtr(p)  = 0._r8
            lgsf(p)  = 0._r8
         end if
      end do
               
   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
   if (CN_evergreen_phenology_opt == 1) then    
   do fp = 1,num_soilp    
      p = filter_soilp(fp)    
      if (evergreen(ivt(p)) == 1._r8) then    
   
         tranr=0.0002_r8   
         ! set carbon fluxes for shifting storage pools to transfer pools    
         leafc_storage_to_xfer(p)  = tranr * leafc_storage(p)/dt    
         frootc_storage_to_xfer(p) = tranr * frootc_storage(p)/dt    
         if (woody(ivt(p)) == 1.0_r8) then    
            livestemc_storage_to_xfer(p)  = tranr * livestemc_storage(p)/dt    
            deadstemc_storage_to_xfer(p)  = tranr * deadstemc_storage(p)/dt    
            livecrootc_storage_to_xfer(p) = tranr * livecrootc_storage(p)/dt   
            deadcrootc_storage_to_xfer(p) = tranr * deadcrootc_storage(p)/dt   
            gresp_storage_to_xfer(p)      = tranr * gresp_storage(p)/dt        
         end if    

        ! set nitrogen fluxes for shifting storage pools to transfer pools    
        leafn_storage_to_xfer(p)  = tranr * leafn_storage(p)/dt    
        frootn_storage_to_xfer(p) = tranr * frootn_storage(p)/dt   
        if (woody(ivt(p)) == 1.0_r8) then    
            livestemn_storage_to_xfer(p)  = tranr * livestemn_storage(p)/dt    
            deadstemn_storage_to_xfer(p)  = tranr * deadstemn_storage(p)/dt    
            livecrootn_storage_to_xfer(p) = tranr * livecrootn_storage(p)/dt   
            deadcrootn_storage_to_xfer(p) = tranr * deadcrootn_storage(p)/dt   
        end if    
                        
        t1 = 1.0_r8 / dt   
            
        leafc_xfer_to_leafc(p)   = t1 * leafc_xfer(p)    
        frootc_xfer_to_frootc(p) = t1 * frootc_xfer(p)   
            
        leafn_xfer_to_leafn(p)   = t1 * leafn_xfer(p)    
        frootn_xfer_to_frootn(p) = t1 * frootn_xfer(p)   
        if (woody(ivt(p)) == 1.0_r8) then   
            livestemc_xfer_to_livestemc(p)   = t1 * livestemc_xfer(p)   
            deadstemc_xfer_to_deadstemc(p)   = t1 * deadstemc_xfer(p)   
            livecrootc_xfer_to_livecrootc(p) = t1 * livecrootc_xfer(p)  
            deadcrootc_xfer_to_deadcrootc(p) = t1 * deadcrootc_xfer(p)  
                
            livestemn_xfer_to_livestemn(p)   = t1 * livestemn_xfer(p)   
            deadstemn_xfer_to_deadstemn(p)   = t1 * deadstemn_xfer(p)   
            livecrootn_xfer_to_livecrootn(p) = t1 * livecrootn_xfer(p)  
            deadcrootn_xfer_to_deadcrootn(p) = t1 * deadcrootn_xfer(p)  
        end if
                
      end if ! end of if (evergreen(ivt(p)) == 1._r8) then    
     
   end do ! end of pft loop 
   
   end if ! end of if (CN_evergreen_phenology_opt == 1) then    
   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    end associate

  end subroutine CNEvergreenPhenology

  !-----------------------------------------------------------------------
  subroutine CNSeasonDecidPhenology (num_soilp, filter_soilp       , &
       temperature_inst, cnveg_state_inst, dgvs_inst , &
       cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! For coupled carbon-nitrogen code (CN).
    ! This routine handles the seasonal deciduous phenology code (temperate
    ! deciduous vegetation that has only one growing season per year).
    !
    ! !USES:
    use shr_const_mod   , only: SHR_CONST_TKFRZ, SHR_CONST_PI
    use clm_varcon      , only: secspday
    use clm_varctl      , only: use_cndv
    !
    ! !ARGUMENTS:
    integer                        , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                        , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(temperature_type)         , intent(in)    :: temperature_inst
    type(cnveg_state_type)         , intent(inout) :: cnveg_state_inst
    type(dgvs_type)                , intent(inout) :: dgvs_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type) , intent(inout) :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    integer :: g,c,p          !indices
    integer :: fp             !lake filter patch index
    real(r8):: ws_flag        !winter-summer solstice flag (0 or 1)
    real(r8):: crit_onset_gdd !critical onset growing degree-day sum
    real(r8):: soilt
    !-----------------------------------------------------------------------

    associate(                                                                                                   & 
         ivt                                 =>    patch%itype                                                   , & ! Input:  [integer   (:)   ]  patch vegetation type                                
         dayl                                =>    grc%dayl                                                    , & ! Input:  [real(r8)  (:)   ]  daylength (s)
         prev_dayl                           =>    grc%prev_dayl                                               , & ! Input:  [real(r8)  (:)   ]  daylength from previous time step (s)
         
         woody                               =>    pftcon%woody                                                , & ! Input:  binary flag for woody lifeform (1=woody, 0=not woody)
         season_decid                        =>    pftcon%season_decid                                         , & ! Input:  binary flag for seasonal-deciduous leaf habit (0 or 1)
         
         t_soisno                            =>    temperature_inst%t_soisno_col                               , & ! Input:  [real(r8)  (:,:) ]  soil temperature (Kelvin)  (-nlevsno+1:nlevgrnd)
         
         pftmayexist                         =>    dgvs_inst%pftmayexist_patch                                 , & ! Output: [logical   (:)   ]  exclude seasonal decid patches from tropics           

         annavg_t2m                          =>    cnveg_state_inst%annavg_t2m_patch                           , & ! Input:  [real(r8)  (:)   ]  annual average 2m air temperature (K)             
         dormant_flag                        =>    cnveg_state_inst%dormant_flag_patch                         , & ! Output: [real(r8)  (:)   ]  dormancy flag                                     
         days_active                         =>    cnveg_state_inst%days_active_patch                          , & ! Output: [real(r8)  (:)   ]  number of days since last dormancy                
         onset_flag                          =>    cnveg_state_inst%onset_flag_patch                           , & ! Output: [real(r8)  (:)   ]  onset flag                                        
         onset_counter                       =>    cnveg_state_inst%onset_counter_patch                        , & ! Output: [real(r8)  (:)   ]  onset counter (seconds)                           
         onset_gddflag                       =>    cnveg_state_inst%onset_gddflag_patch                        , & ! Output: [real(r8)  (:)   ]  onset freeze flag                                 
         onset_gdd                           =>    cnveg_state_inst%onset_gdd_patch                            , & ! Output: [real(r8)  (:)   ]  onset growing degree days                         
         offset_flag                         =>    cnveg_state_inst%offset_flag_patch                          , & ! Output: [real(r8)  (:)   ]  offset flag                                       
         offset_counter                      =>    cnveg_state_inst%offset_counter_patch                       , & ! Output: [real(r8)  (:)   ]  offset counter (seconds)                          
         bglfr                               =>    cnveg_state_inst%bglfr_patch                                , & ! Output: [real(r8)  (:)   ]  background litterfall rate (1/s)                  
         bgtr                                =>    cnveg_state_inst%bgtr_patch                                 , & ! Output: [real(r8)  (:)   ]  background transfer growth rate (1/s)             
         lgsf                                =>    cnveg_state_inst%lgsf_patch                                 , & ! Output: [real(r8)  (:)   ]  long growing season factor [0-1]                  
         
         leafc_storage                       =>    cnveg_carbonstate_inst%leafc_storage_patch                  , & ! Input:  [real(r8)  (:)   ]  (gC/m2) leaf C storage                            
         frootc_storage                      =>    cnveg_carbonstate_inst%frootc_storage_patch                 , & ! Input:  [real(r8)  (:)   ]  (gC/m2) fine root C storage                       
         livestemc_storage                   =>    cnveg_carbonstate_inst%livestemc_storage_patch              , & ! Input:  [real(r8)  (:)   ]  (gC/m2) live stem C storage                       
         deadstemc_storage                   =>    cnveg_carbonstate_inst%deadstemc_storage_patch              , & ! Input:  [real(r8)  (:)   ]  (gC/m2) dead stem C storage                       
         livecrootc_storage                  =>    cnveg_carbonstate_inst%livecrootc_storage_patch             , & ! Input:  [real(r8)  (:)   ]  (gC/m2) live coarse root C storage                
         deadcrootc_storage                  =>    cnveg_carbonstate_inst%deadcrootc_storage_patch             , & ! Input:  [real(r8)  (:)   ]  (gC/m2) dead coarse root C storage                
         gresp_storage                       =>    cnveg_carbonstate_inst%gresp_storage_patch                  , & ! Input:  [real(r8)  (:)   ]  (gC/m2) growth respiration storage                
         leafc_xfer                          =>    cnveg_carbonstate_inst%leafc_xfer_patch                     , & ! Output:  [real(r8) (:)   ]  (gC/m2) leaf C transfer                           
         frootc_xfer                         =>    cnveg_carbonstate_inst%frootc_xfer_patch                    , & ! Output:  [real(r8) (:)   ]  (gC/m2) fine root C transfer                      
         livestemc_xfer                      =>    cnveg_carbonstate_inst%livestemc_xfer_patch                 , & ! Output:  [real(r8) (:)   ]  (gC/m2) live stem C transfer                      
         deadstemc_xfer                      =>    cnveg_carbonstate_inst%deadstemc_xfer_patch                 , & ! Output:  [real(r8) (:)   ]  (gC/m2) dead stem C transfer                      
         livecrootc_xfer                     =>    cnveg_carbonstate_inst%livecrootc_xfer_patch                , & ! Output:  [real(r8) (:)   ]  (gC/m2) live coarse root C transfer               
         deadcrootc_xfer                     =>    cnveg_carbonstate_inst%deadcrootc_xfer_patch                , & ! Output:  [real(r8) (:)   ]  (gC/m2) dead coarse root C transfer               
         
         leafn_storage                       =>    cnveg_nitrogenstate_inst%leafn_storage_patch                , & ! Input:  [real(r8)  (:)   ]  (gN/m2) leaf N storage                            
         frootn_storage                      =>    cnveg_nitrogenstate_inst%frootn_storage_patch               , & ! Input:  [real(r8)  (:)   ]  (gN/m2) fine root N storage                       
         livestemn_storage                   =>    cnveg_nitrogenstate_inst%livestemn_storage_patch            , & ! Input:  [real(r8)  (:)   ]  (gN/m2) live stem N storage                       
         deadstemn_storage                   =>    cnveg_nitrogenstate_inst%deadstemn_storage_patch            , & ! Input:  [real(r8)  (:)   ]  (gN/m2) dead stem N storage                       
         livecrootn_storage                  =>    cnveg_nitrogenstate_inst%livecrootn_storage_patch           , & ! Input:  [real(r8)  (:)   ]  (gN/m2) live coarse root N storage                
         deadcrootn_storage                  =>    cnveg_nitrogenstate_inst%deadcrootn_storage_patch           , & ! Input:  [real(r8)  (:)   ]  (gN/m2) dead coarse root N storage                
         leafn_xfer                          =>    cnveg_nitrogenstate_inst%leafn_xfer_patch                   , & ! Output:  [real(r8) (:)   ]  (gN/m2) leaf N transfer                           
         frootn_xfer                         =>    cnveg_nitrogenstate_inst%frootn_xfer_patch                  , & ! Output:  [real(r8) (:)   ]  (gN/m2) fine root N transfer                      
         livestemn_xfer                      =>    cnveg_nitrogenstate_inst%livestemn_xfer_patch               , & ! Output:  [real(r8) (:)   ]  (gN/m2) live stem N transfer                      
         deadstemn_xfer                      =>    cnveg_nitrogenstate_inst%deadstemn_xfer_patch               , & ! Output:  [real(r8) (:)   ]  (gN/m2) dead stem N transfer                      
         livecrootn_xfer                     =>    cnveg_nitrogenstate_inst%livecrootn_xfer_patch              , & ! Output:  [real(r8) (:)   ]  (gN/m2) live coarse root N transfer               
         deadcrootn_xfer                     =>    cnveg_nitrogenstate_inst%deadcrootn_xfer_patch              , & ! Output:  [real(r8) (:)   ]  (gN/m2) dead coarse root N transfer               

         prev_leafc_to_litter                =>    cnveg_carbonflux_inst%prev_leafc_to_litter_patch            , & ! Output: [real(r8)  (:)   ]  previous timestep leaf C litterfall flux (gC/m2/s)
         prev_frootc_to_litter               =>    cnveg_carbonflux_inst%prev_frootc_to_litter_patch           , & ! Output: [real(r8)  (:)   ]  previous timestep froot C litterfall flux (gC/m2/s)
         leafc_xfer_to_leafc                 =>    cnveg_carbonflux_inst%leafc_xfer_to_leafc_patch             , & ! Output:  [real(r8) (:)   ]                                                    
         frootc_xfer_to_frootc               =>    cnveg_carbonflux_inst%frootc_xfer_to_frootc_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         livestemc_xfer_to_livestemc         =>    cnveg_carbonflux_inst%livestemc_xfer_to_livestemc_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemc_xfer_to_deadstemc         =>    cnveg_carbonflux_inst%deadstemc_xfer_to_deadstemc_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootc_xfer_to_livecrootc       =>    cnveg_carbonflux_inst%livecrootc_xfer_to_livecrootc_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootc_xfer_to_deadcrootc       =>    cnveg_carbonflux_inst%deadcrootc_xfer_to_deadcrootc_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         leafc_storage_to_xfer               =>    cnveg_carbonflux_inst%leafc_storage_to_xfer_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         frootc_storage_to_xfer              =>    cnveg_carbonflux_inst%frootc_storage_to_xfer_patch          , & ! Output:  [real(r8) (:)   ]                                                    
         livestemc_storage_to_xfer           =>    cnveg_carbonflux_inst%livestemc_storage_to_xfer_patch       , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemc_storage_to_xfer           =>    cnveg_carbonflux_inst%deadstemc_storage_to_xfer_patch       , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootc_storage_to_xfer          =>    cnveg_carbonflux_inst%livecrootc_storage_to_xfer_patch      , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootc_storage_to_xfer          =>    cnveg_carbonflux_inst%deadcrootc_storage_to_xfer_patch      , & ! Output:  [real(r8) (:)   ]                                                    
         gresp_storage_to_xfer               =>    cnveg_carbonflux_inst%gresp_storage_to_xfer_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         
         leafn_xfer_to_leafn                 =>    cnveg_nitrogenflux_inst%leafn_xfer_to_leafn_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         frootn_xfer_to_frootn               =>    cnveg_nitrogenflux_inst%frootn_xfer_to_frootn_patch         , & ! Output:  [real(r8) (:)   ]                                                    
         livestemn_xfer_to_livestemn         =>    cnveg_nitrogenflux_inst%livestemn_xfer_to_livestemn_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemn_xfer_to_deadstemn         =>    cnveg_nitrogenflux_inst%deadstemn_xfer_to_deadstemn_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootn_xfer_to_livecrootn       =>    cnveg_nitrogenflux_inst%livecrootn_xfer_to_livecrootn_patch , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootn_xfer_to_deadcrootn       =>    cnveg_nitrogenflux_inst%deadcrootn_xfer_to_deadcrootn_patch , & ! Output:  [real(r8) (:)   ]                                                    
         leafn_storage_to_xfer               =>    cnveg_nitrogenflux_inst%leafn_storage_to_xfer_patch         , & ! Output:  [real(r8) (:)   ]                                                    
         frootn_storage_to_xfer              =>    cnveg_nitrogenflux_inst%frootn_storage_to_xfer_patch        , & ! Output:  [real(r8) (:)   ]                                                    
         livestemn_storage_to_xfer           =>    cnveg_nitrogenflux_inst%livestemn_storage_to_xfer_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemn_storage_to_xfer           =>    cnveg_nitrogenflux_inst%deadstemn_storage_to_xfer_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootn_storage_to_xfer          =>    cnveg_nitrogenflux_inst%livecrootn_storage_to_xfer_patch    , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootn_storage_to_xfer          =>    cnveg_nitrogenflux_inst%deadcrootn_storage_to_xfer_patch      & ! Output:  [real(r8) (:)   ]                                                    
         )

      ! start patch loop
      do fp = 1,num_soilp
         p = filter_soilp(fp)
         c = patch%column(p)
         g = patch%gridcell(p)

         if (season_decid(ivt(p)) == 1._r8) then

            ! set background litterfall rate, background transfer rate, and
            ! long growing season factor to 0 for seasonal deciduous types
            bglfr(p) = 0._r8
            bgtr(p) = 0._r8
            lgsf(p) = 0._r8

            ! onset gdd sum from Biome-BGC, v4.1.2
            crit_onset_gdd = exp(4.8_r8 + 0.13_r8*(annavg_t2m(p) - SHR_CONST_TKFRZ))

            ! set flag for solstice period (winter->summer = 1, summer->winter = 0)
            if (dayl(g) >= prev_dayl(g)) then
               ws_flag = 1._r8
            else
               ws_flag = 0._r8
            end if

            ! update offset_counter and test for the end of the offset period
            if (offset_flag(p) == 1.0_r8) then
               ! decrement counter for offset period
               offset_counter(p) = offset_counter(p) - dt

               ! if this is the end of the offset_period, reset phenology
               ! flags and indices
               if (offset_counter(p) == 0.0_r8) then
                  ! this code block was originally handled by call cn_offset_cleanup(p)
                  ! inlined during vectorization

                  offset_flag(p) = 0._r8
                  offset_counter(p) = 0._r8
                  dormant_flag(p) = 1._r8
                  days_active(p) = 0._r8
                  if (use_cndv) then
                     pftmayexist(p) = .true.
                  end if

                  ! reset the previous timestep litterfall flux memory
                  prev_leafc_to_litter(p) = 0._r8
                  prev_frootc_to_litter(p) = 0._r8
               end if
            end if

            ! update onset_counter and test for the end of the onset period
            if (onset_flag(p) == 1.0_r8) then
               ! decrement counter for onset period
               onset_counter(p) = onset_counter(p) - dt

               ! if this is the end of the onset period, reset phenology
               ! flags and indices
               if (onset_counter(p) == 0.0_r8) then
                  ! this code block was originally handled by call cn_onset_cleanup(p)
                  ! inlined during vectorization

                  onset_flag(p) = 0.0_r8
                  onset_counter(p) = 0.0_r8
                  ! set all transfer growth rates to 0.0
                  leafc_xfer_to_leafc(p)   = 0.0_r8
                  frootc_xfer_to_frootc(p) = 0.0_r8
                  leafn_xfer_to_leafn(p)   = 0.0_r8
                  frootn_xfer_to_frootn(p) = 0.0_r8
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemc_xfer_to_livestemc(p)   = 0.0_r8
                     deadstemc_xfer_to_deadstemc(p)   = 0.0_r8
                     livecrootc_xfer_to_livecrootc(p) = 0.0_r8
                     deadcrootc_xfer_to_deadcrootc(p) = 0.0_r8
                     livestemn_xfer_to_livestemn(p)   = 0.0_r8
                     deadstemn_xfer_to_deadstemn(p)   = 0.0_r8
                     livecrootn_xfer_to_livecrootn(p) = 0.0_r8
                     deadcrootn_xfer_to_deadcrootn(p) = 0.0_r8
                  end if
                  ! set transfer pools to 0.0
                  leafc_xfer(p) = 0.0_r8
                  leafn_xfer(p) = 0.0_r8
                  frootc_xfer(p) = 0.0_r8
                  frootn_xfer(p) = 0.0_r8
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemc_xfer(p) = 0.0_r8
                     livestemn_xfer(p) = 0.0_r8
                     deadstemc_xfer(p) = 0.0_r8
                     deadstemn_xfer(p) = 0.0_r8
                     livecrootc_xfer(p) = 0.0_r8
                     livecrootn_xfer(p) = 0.0_r8
                     deadcrootc_xfer(p) = 0.0_r8
                     deadcrootn_xfer(p) = 0.0_r8
                  end if
               end if
            end if

            ! test for switching from dormant period to growth period
            if (dormant_flag(p) == 1.0_r8) then

               ! Test to turn on growing degree-day sum, if off.
               ! switch on the growing degree day sum on the winter solstice

               if (onset_gddflag(p) == 0._r8 .and. ws_flag == 1._r8) then
                  onset_gddflag(p) = 1._r8
                  onset_gdd(p) = 0._r8
               end if

               ! Test to turn off growing degree-day sum, if on.
               ! This test resets the growing degree day sum if it gets past
               ! the summer solstice without reaching the threshold value.
               ! In that case, it will take until the next winter solstice
               ! before the growing degree-day summation starts again.

               if (onset_gddflag(p) == 1._r8 .and. ws_flag == 0._r8) then
                  onset_gddflag(p) = 0._r8
                  onset_gdd(p) = 0._r8
               end if

               ! if the gdd flag is set, and if the soil is above freezing
               ! then accumulate growing degree days for onset trigger

               soilt = t_soisno(c, phenology_soil_layer)
               if (onset_gddflag(p) == 1.0_r8 .and. soilt > SHR_CONST_TKFRZ) then
                  onset_gdd(p) = onset_gdd(p) + (soilt-SHR_CONST_TKFRZ)*fracday
               end if

               ! set onset_flag if critical growing degree-day sum is exceeded
               if (onset_gdd(p) > crit_onset_gdd) then
                  onset_flag(p) = 1.0_r8
                  dormant_flag(p) = 0.0_r8
                  onset_gddflag(p) = 0.0_r8
                  onset_gdd(p) = 0.0_r8
                  onset_counter(p) = ndays_on * secspday

                  ! move all the storage pools into transfer pools,
                  ! where they will be transfered to displayed growth over the onset period.
                  ! this code was originally handled with call cn_storage_to_xfer(p)
                  ! inlined during vectorization

                  ! set carbon fluxes for shifting storage pools to transfer pools
                  leafc_storage_to_xfer(p)  = fstor2tran * leafc_storage(p)/dt
                  frootc_storage_to_xfer(p) = fstor2tran * frootc_storage(p)/dt
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemc_storage_to_xfer(p)  = fstor2tran * livestemc_storage(p)/dt
                     deadstemc_storage_to_xfer(p)  = fstor2tran * deadstemc_storage(p)/dt
                     livecrootc_storage_to_xfer(p) = fstor2tran * livecrootc_storage(p)/dt
                     deadcrootc_storage_to_xfer(p) = fstor2tran * deadcrootc_storage(p)/dt
                     gresp_storage_to_xfer(p)      = fstor2tran * gresp_storage(p)/dt
                  end if

                  ! set nitrogen fluxes for shifting storage pools to transfer pools
                  leafn_storage_to_xfer(p)  = fstor2tran * leafn_storage(p)/dt
                  frootn_storage_to_xfer(p) = fstor2tran * frootn_storage(p)/dt
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemn_storage_to_xfer(p)  = fstor2tran * livestemn_storage(p)/dt
                     deadstemn_storage_to_xfer(p)  = fstor2tran * deadstemn_storage(p)/dt
                     livecrootn_storage_to_xfer(p) = fstor2tran * livecrootn_storage(p)/dt
                     deadcrootn_storage_to_xfer(p) = fstor2tran * deadcrootn_storage(p)/dt
                  end if
               end if

               ! test for switching from growth period to offset period
            else if (offset_flag(p) == 0.0_r8) then
               if (use_cndv) then
                  ! If days_active > 355, then remove patch in
                  ! CNDVEstablishment at the end of the year.
                  ! days_active > 355 is a symptom of seasonal decid. patches occurring in
                  ! gridcells where dayl never drops below crit_dayl.
                  ! This results in TLAI>1e4 in a few gridcells.
                  days_active(p) = days_active(p) + fracday
                  if (days_active(p) > 355._r8) pftmayexist(p) = .false.
               end if

               ! only begin to test for offset daylength once past the summer sol
               if (ws_flag == 0._r8 .and. dayl(g) < crit_dayl) then
                  offset_flag(p) = 1._r8
                  offset_counter(p) = ndays_off * secspday
                  prev_leafc_to_litter(p) = 0._r8
                  prev_frootc_to_litter(p) = 0._r8
               end if
            end if

         end if ! end if seasonal deciduous

      end do ! end of patch loop

    end associate
 
  end subroutine CNSeasonDecidPhenology

  !-----------------------------------------------------------------------
  subroutine CNStressDecidPhenology (num_soilp, filter_soilp , &
       soilstate_inst, temperature_inst, atm2lnd_inst, wateratm2lndbulk_inst, cnveg_state_inst, &
       cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, &
       cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! This routine handles phenology for vegetation types, such as grasses and
    ! tropical drought deciduous trees, that respond to cold and drought stress
    ! signals and that can have multiple growing seasons in a given year.
    ! This routine allows for the possibility that leaves might persist year-round
    ! in the absence of a suitable stress trigger, by switching to an essentially
    ! evergreen habit, but maintaining a deciduous leaf longevity, while waiting
    ! for the next stress trigger.  This is in contrast to the seasonal deciduous
    ! algorithm (for temperate deciduous trees) that forces a single growing season
    ! per year.
    !
    ! !USES:
    use clm_time_manager , only : get_days_per_year
    use CNSharedParamsMod, only : use_fun
    use clm_varcon       , only : secspday
    use shr_const_mod    , only : SHR_CONST_TKFRZ, SHR_CONST_PI
    use CNSharedParamsMod, only : CNParamsShareInst
    !
    ! !ARGUMENTS:
    integer                        , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                        , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(soilstate_type)           , intent(in)    :: soilstate_inst
    type(temperature_type)         , intent(in)    :: temperature_inst
    type(atm2lnd_type)             , intent(in)    :: atm2lnd_inst
    type(wateratm2lndbulk_type)             , intent(in)    :: wateratm2lndbulk_inst
    type(cnveg_state_type)         , intent(inout) :: cnveg_state_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type) , intent(inout) :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    real(r8),parameter :: secspqtrday = secspday / 4  ! seconds per quarter day
    integer :: g,c,p           ! indices
    integer :: fp              ! lake filter patch index
    real(r8):: dayspyr         ! days per year
    real(r8):: crit_onset_gdd  ! degree days for onset trigger
    real(r8):: soilt           ! temperature of top soil layer
    real(r8):: psi             ! water stress of top soil layer
    real(r8):: rain_threshold  ! rain threshold for leaf on [mm]
    logical :: additional_onset_condition ! additional condition for leaf onset
    !-----------------------------------------------------------------------

    associate(                                                                                                   & 
         ivt                                 =>    patch%itype                                                   , & ! Input:  [integer   (:)   ]  patch vegetation type                                
         dayl                                =>    grc%dayl                                                    , & ! Input:  [real(r8)  (:)   ]  daylength (s)
         
         prec10                              => wateratm2lndbulk_inst%prec10_patch                                     , & ! Input:  [real(r8) (:)     ]  10-day running mean of tot. precipitation
         leaf_long                           =>    pftcon%leaf_long                                            , & ! Input:  leaf longevity (yrs)                              
         woody                               =>    pftcon%woody                                                , & ! Input:  binary flag for woody lifeform (1=woody, 0=not woody)
         stress_decid                        =>    pftcon%stress_decid                                         , & ! Input:  binary flag for stress-deciduous leaf habit (0 or 1)
         
         soilpsi                             =>    soilstate_inst%soilpsi_col                                  , & ! Input:  [real(r8)  (:,:) ]  soil water potential in each soil layer (MPa)   
         
         t_soisno                            =>    temperature_inst%t_soisno_col                               , & ! Input:  [real(r8)  (:,:) ]  soil temperature (Kelvin)  (-nlevsno+1:nlevgrnd)
         dormant_flag                        =>    cnveg_state_inst%dormant_flag_patch                         , & ! Output:  [real(r8) (:)   ]  dormancy flag                                     
         days_active                         =>    cnveg_state_inst%days_active_patch                          , & ! Output:  [real(r8) (:)   ]  number of days since last dormancy                
         onset_flag                          =>    cnveg_state_inst%onset_flag_patch                           , & ! Output:  [real(r8) (:)   ]  onset flag                                        
         onset_counter                       =>    cnveg_state_inst%onset_counter_patch                        , & ! Output:  [real(r8) (:)   ]  onset counter (seconds)                           
         onset_gddflag                       =>    cnveg_state_inst%onset_gddflag_patch                        , & ! Output:  [real(r8) (:)   ]  onset freeze flag                                 
         onset_fdd                           =>    cnveg_state_inst%onset_fdd_patch                            , & ! Output:  [real(r8) (:)   ]  onset freezing degree days counter                
         onset_gdd                           =>    cnveg_state_inst%onset_gdd_patch                            , & ! Output:  [real(r8) (:)   ]  onset growing degree days                         
         onset_swi                           =>    cnveg_state_inst%onset_swi_patch                            , & ! Output:  [real(r8) (:)   ]  onset soil water index                            
         offset_flag                         =>    cnveg_state_inst%offset_flag_patch                          , & ! Output:  [real(r8) (:)   ]  offset flag                                       
         offset_counter                      =>    cnveg_state_inst%offset_counter_patch                       , & ! Output:  [real(r8) (:)   ]  offset counter (seconds)                          
         offset_fdd                          =>    cnveg_state_inst%offset_fdd_patch                           , & ! Output:  [real(r8) (:)   ]  offset freezing degree days counter               
         offset_swi                          =>    cnveg_state_inst%offset_swi_patch                           , & ! Output:  [real(r8) (:)   ]  offset soil water index                           
         lgsf                                =>    cnveg_state_inst%lgsf_patch                                 , & ! Output:  [real(r8) (:)   ]  long growing season factor [0-1]                  
         bglfr                               =>    cnveg_state_inst%bglfr_patch                                , & ! Output:  [real(r8) (:)   ]  background litterfall rate (1/s)                  
         bgtr                                =>    cnveg_state_inst%bgtr_patch                                 , & ! Output:  [real(r8) (:)   ]  background transfer growth rate (1/s)             
         annavg_t2m                          =>    cnveg_state_inst%annavg_t2m_patch                           , & ! Output:  [real(r8) (:)   ]  annual average 2m air temperature (K)             
         leafc                               =>    cnveg_carbonstate_inst%leafc_patch                          , & ! Input:  [real(r8)  (:)   ]  (gC/m2) leaf C 
     
         frootc                              =>    cnveg_carbonstate_inst%frootc_patch                         , & ! Input:  [real(r8) (:)    ]  (gC/m2) fine root C
         leafc_storage                       =>    cnveg_carbonstate_inst%leafc_storage_patch                  , & ! Input:  [real(r8)  (:)   ]  (gC/m2) leaf C storage                            
         frootc_storage                      =>    cnveg_carbonstate_inst%frootc_storage_patch                 , & ! Input:  [real(r8)  (:)   ]  (gC/m2) fine root C storage                       
         livestemc_storage                   =>    cnveg_carbonstate_inst%livestemc_storage_patch              , & ! Input:  [real(r8)  (:)   ]  (gC/m2) live stem C storage                       
         deadstemc_storage                   =>    cnveg_carbonstate_inst%deadstemc_storage_patch              , & ! Input:  [real(r8)  (:)   ]  (gC/m2) dead stem C storage                       
         livecrootc_storage                  =>    cnveg_carbonstate_inst%livecrootc_storage_patch             , & ! Input:  [real(r8)  (:)   ]  (gC/m2) live coarse root C storage                
         deadcrootc_storage                  =>    cnveg_carbonstate_inst%deadcrootc_storage_patch             , & ! Input:  [real(r8)  (:)   ]  (gC/m2) dead coarse root C storage                
         gresp_storage                       =>    cnveg_carbonstate_inst%gresp_storage_patch                  , & ! Input:  [real(r8)  (:)   ]  (gC/m2) growth respiration storage                
         leafc_xfer                          =>    cnveg_carbonstate_inst%leafc_xfer_patch                     , & ! Output:  [real(r8) (:)   ]  (gC/m2) leaf C transfer
         frootc_xfer                         =>    cnveg_carbonstate_inst%frootc_xfer_patch                    , & ! Output:  [real(r8) (:)   ]  (gC/m2) fine root C transfer                      
         livestemc_xfer                      =>    cnveg_carbonstate_inst%livestemc_xfer_patch                 , & ! Output:  [real(r8) (:)   ]  (gC/m2) live stem C transfer                      
         deadstemc_xfer                      =>    cnveg_carbonstate_inst%deadstemc_xfer_patch                 , & ! Output:  [real(r8) (:)   ]  (gC/m2) dead stem C transfer                      
         livecrootc_xfer                     =>    cnveg_carbonstate_inst%livecrootc_xfer_patch                , & ! Output:  [real(r8) (:)   ]  (gC/m2) live coarse root C transfer               
         deadcrootc_xfer                     =>    cnveg_carbonstate_inst%deadcrootc_xfer_patch                , & ! Output:  [real(r8) (:)   ]  (gC/m2) dead coarse root C transfer               
         leafn_storage                       =>    cnveg_nitrogenstate_inst%leafn_storage_patch                , & ! Input:  [real(r8)  (:)   ]  (gN/m2) leaf N storage                            
         frootn_storage                      =>    cnveg_nitrogenstate_inst%frootn_storage_patch               , & ! Input:  [real(r8)  (:)   ]  (gN/m2) fine root N storage                       
         livestemn_storage                   =>    cnveg_nitrogenstate_inst%livestemn_storage_patch            , & ! Input:  [real(r8)  (:)   ]  (gN/m2) live stem N storage                       
         deadstemn_storage                   =>    cnveg_nitrogenstate_inst%deadstemn_storage_patch            , & ! Input:  [real(r8)  (:)   ]  (gN/m2) dead stem N storage                       
         livecrootn_storage                  =>    cnveg_nitrogenstate_inst%livecrootn_storage_patch           , & ! Input:  [real(r8)  (:)   ]  (gN/m2) live coarse root N storage                
         deadcrootn_storage                  =>    cnveg_nitrogenstate_inst%deadcrootn_storage_patch           , & ! Input:  [real(r8)  (:)   ]  (gN/m2) dead coarse root N storage                
         leafn_xfer                          =>    cnveg_nitrogenstate_inst%leafn_xfer_patch                   , & ! Output:  [real(r8) (:)   ]  (gN/m2) leaf N transfer                           
         frootn_xfer                         =>    cnveg_nitrogenstate_inst%frootn_xfer_patch                  , & ! Output:  [real(r8) (:)   ]  (gN/m2) fine root N transfer                      
         livestemn_xfer                      =>    cnveg_nitrogenstate_inst%livestemn_xfer_patch               , & ! Output:  [real(r8) (:)   ]  (gN/m2) live stem N transfer                      
         deadstemn_xfer                      =>    cnveg_nitrogenstate_inst%deadstemn_xfer_patch               , & ! Output:  [real(r8) (:)   ]  (gN/m2) dead stem N transfer                      
         livecrootn_xfer                     =>    cnveg_nitrogenstate_inst%livecrootn_xfer_patch              , & ! Output:  [real(r8) (:)   ]  (gN/m2) live coarse root N transfer               
         deadcrootn_xfer                     =>    cnveg_nitrogenstate_inst%deadcrootn_xfer_patch              , & ! Output:  [real(r8) (:)   ]  (gN/m2) dead coarse root N transfer               
         
         prev_leafc_to_litter                =>    cnveg_carbonflux_inst%prev_leafc_to_litter_patch            , & ! Output:  [real(r8) (:)   ]  previous timestep leaf C litterfall flux (gC/m2/s)
         prev_frootc_to_litter               =>    cnveg_carbonflux_inst%prev_frootc_to_litter_patch           , & ! Output:  [real(r8) (:)   ]  previous timestep froot C litterfall flux (gC/m2/s)
         leafc_xfer_to_leafc                 =>    cnveg_carbonflux_inst%leafc_xfer_to_leafc_patch             , & ! Output:  [real(r8) (:)   ]                                                    
         frootc_xfer_to_frootc               =>    cnveg_carbonflux_inst%frootc_xfer_to_frootc_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         livestemc_xfer_to_livestemc         =>    cnveg_carbonflux_inst%livestemc_xfer_to_livestemc_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemc_xfer_to_deadstemc         =>    cnveg_carbonflux_inst%deadstemc_xfer_to_deadstemc_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootc_xfer_to_livecrootc       =>    cnveg_carbonflux_inst%livecrootc_xfer_to_livecrootc_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootc_xfer_to_deadcrootc       =>    cnveg_carbonflux_inst%deadcrootc_xfer_to_deadcrootc_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         leafc_storage_to_xfer               =>    cnveg_carbonflux_inst%leafc_storage_to_xfer_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         frootc_storage_to_xfer              =>    cnveg_carbonflux_inst%frootc_storage_to_xfer_patch          , & ! Output:  [real(r8) (:)   ]                                                    
         livestemc_storage_to_xfer           =>    cnveg_carbonflux_inst%livestemc_storage_to_xfer_patch       , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemc_storage_to_xfer           =>    cnveg_carbonflux_inst%deadstemc_storage_to_xfer_patch       , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootc_storage_to_xfer          =>    cnveg_carbonflux_inst%livecrootc_storage_to_xfer_patch      , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootc_storage_to_xfer          =>    cnveg_carbonflux_inst%deadcrootc_storage_to_xfer_patch      , & ! Output:  [real(r8) (:)   ]                                                    
         gresp_storage_to_xfer               =>    cnveg_carbonflux_inst%gresp_storage_to_xfer_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         
         leafn_xfer_to_leafn                 =>    cnveg_nitrogenflux_inst%leafn_xfer_to_leafn_patch           , & ! Output:  [real(r8) (:)   ]                                                    
         frootn_xfer_to_frootn               =>    cnveg_nitrogenflux_inst%frootn_xfer_to_frootn_patch         , & ! Output:  [real(r8) (:)   ]                                                    
         livestemn_xfer_to_livestemn         =>    cnveg_nitrogenflux_inst%livestemn_xfer_to_livestemn_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemn_xfer_to_deadstemn         =>    cnveg_nitrogenflux_inst%deadstemn_xfer_to_deadstemn_patch   , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootn_xfer_to_livecrootn       =>    cnveg_nitrogenflux_inst%livecrootn_xfer_to_livecrootn_patch , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootn_xfer_to_deadcrootn       =>    cnveg_nitrogenflux_inst%deadcrootn_xfer_to_deadcrootn_patch , & ! Output:  [real(r8) (:)   ]                                                    
         leafn_storage_to_xfer               =>    cnveg_nitrogenflux_inst%leafn_storage_to_xfer_patch         , & ! Output:  [real(r8) (:)   ]                                                    
         frootn_storage_to_xfer              =>    cnveg_nitrogenflux_inst%frootn_storage_to_xfer_patch        , & ! Output:  [real(r8) (:)   ]                                                    
         livestemn_storage_to_xfer           =>    cnveg_nitrogenflux_inst%livestemn_storage_to_xfer_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         deadstemn_storage_to_xfer           =>    cnveg_nitrogenflux_inst%deadstemn_storage_to_xfer_patch     , & ! Output:  [real(r8) (:)   ]                                                    
         livecrootn_storage_to_xfer          =>    cnveg_nitrogenflux_inst%livecrootn_storage_to_xfer_patch    , & ! Output:  [real(r8) (:)   ]                                                    
         deadcrootn_storage_to_xfer          =>    cnveg_nitrogenflux_inst%deadcrootn_storage_to_xfer_patch      & ! Output:  [real(r8) (:)   ]                                                    
         )

      ! set time steps
      dayspyr = get_days_per_year()

      ! specify rain threshold for leaf onset
      rain_threshold = 20._r8

      do fp = 1,num_soilp
         p = filter_soilp(fp)
         c = patch%column(p)
         g = patch%gridcell(p)

         if (stress_decid(ivt(p)) == 1._r8) then
            soilt = t_soisno(c, phenology_soil_layer)
            psi = soilpsi(c, phenology_soil_layer)

            ! onset gdd sum from Biome-BGC, v4.1.2
            crit_onset_gdd = exp(4.8_r8 + 0.13_r8*(annavg_t2m(p) - SHR_CONST_TKFRZ))


            ! update offset_counter and test for the end of the offset period
            if (offset_flag(p) == 1._r8) then
               ! decrement counter for offset period
               offset_counter(p) = offset_counter(p) - dt

               ! if this is the end of the offset_period, reset phenology
               ! flags and indices
               if (offset_counter(p) == 0._r8) then
                  ! this code block was originally handled by call cn_offset_cleanup(p)
                  ! inlined during vectorization
                  offset_flag(p) = 0._r8
                  offset_counter(p) = 0._r8
                  dormant_flag(p) = 1._r8
                  days_active(p) = 0._r8

                  ! reset the previous timestep litterfall flux memory
                  prev_leafc_to_litter(p) = 0._r8
                  prev_frootc_to_litter(p) = 0._r8
               end if
            end if

            ! update onset_counter and test for the end of the onset period
            if (onset_flag(p) == 1.0_r8) then
               ! decrement counter for onset period
               onset_counter(p) = onset_counter(p) - dt

               ! if this is the end of the onset period, reset phenology
               ! flags and indices
               if (onset_counter(p) == 0.0_r8) then
                  ! this code block was originally handled by call cn_onset_cleanup(p)
                  ! inlined during vectorization
                  onset_flag(p) = 0._r8
                  onset_counter(p) = 0._r8
                  ! set all transfer growth rates to 0.0
                  leafc_xfer_to_leafc(p)   = 0._r8
                  frootc_xfer_to_frootc(p) = 0._r8
                  leafn_xfer_to_leafn(p)   = 0._r8
                  frootn_xfer_to_frootn(p) = 0._r8
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemc_xfer_to_livestemc(p)   = 0._r8
                     deadstemc_xfer_to_deadstemc(p)   = 0._r8
                     livecrootc_xfer_to_livecrootc(p) = 0._r8
                     deadcrootc_xfer_to_deadcrootc(p) = 0._r8
                     livestemn_xfer_to_livestemn(p)   = 0._r8
                     deadstemn_xfer_to_deadstemn(p)   = 0._r8
                     livecrootn_xfer_to_livecrootn(p) = 0._r8
                     deadcrootn_xfer_to_deadcrootn(p) = 0._r8
                  end if
                  ! set transfer pools to 0.0
                  leafc_xfer(p) = 0._r8
                  leafn_xfer(p) = 0._r8
                  frootc_xfer(p) = 0._r8
                  frootn_xfer(p) = 0._r8
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemc_xfer(p) = 0._r8
                     livestemn_xfer(p) = 0._r8
                     deadstemc_xfer(p) = 0._r8
                     deadstemn_xfer(p) = 0._r8
                     livecrootc_xfer(p) = 0._r8
                     livecrootn_xfer(p) = 0._r8
                     deadcrootc_xfer(p) = 0._r8
                     deadcrootn_xfer(p) = 0._r8
                  end if
               end if
            end if

            ! test for switching from dormant period to growth period
            if (dormant_flag(p) == 1._r8) then

               ! keep track of the number of freezing degree days in this
               ! dormancy period (only if the freeze flag has not previously been set
               ! for this dormancy period

               if (onset_gddflag(p) == 0._r8 .and. soilt < SHR_CONST_TKFRZ) onset_fdd(p) = onset_fdd(p) + fracday

               ! if the number of freezing degree days exceeds a critical value,
               ! then onset will require both wet soils and a critical soil
               ! temperature sum.  If this case is triggered, reset any previously
               ! accumulated value in onset_swi, so that onset now depends on
               ! the accumulated soil water index following the freeze trigger

               if (onset_fdd(p) > crit_onset_fdd) then
                  onset_gddflag(p) = 1._r8
                  onset_fdd(p) = 0._r8
                  onset_swi(p) = 0._r8
               end if

               ! if the freeze flag is set, and if the soil is above freezing
               ! then accumulate growing degree days for onset trigger

               if (onset_gddflag(p) == 1._r8 .and. soilt > SHR_CONST_TKFRZ) then
                  onset_gdd(p) = onset_gdd(p) + (soilt-SHR_CONST_TKFRZ)*fracday
               end if

              ! if soils are wet, accumulate soil water index for onset trigger
               additional_onset_condition = .true.
               if(CNParamsShareInst%constrain_stress_deciduous_onset) then
                  ! if additional constraint condition not met,  set to false
                  if ((prec10(p) * (3600.0_r8*10.0_r8*24.0_r8)) < rain_threshold) then
                     additional_onset_condition = .false.
                  endif
               endif

               if (psi >= soilpsi_on) then
                  onset_swi(p) = onset_swi(p) + fracday
               endif

               ! if critical soil water index is exceeded, set onset_flag, and
               ! then test for soil temperature criteria

               ! Adding in Kyla's rainfall trigger when fun on. RF. prec10 (mm/s) needs to be higher than 8mm over 10 days. 

               if (onset_swi(p) > crit_onset_swi.and. additional_onset_condition)  then
                  onset_flag(p) = 1._r8
              
                  ! only check soil temperature criteria if freeze flag set since
                  ! beginning of last dormancy.  If freeze flag set and growing
                  ! degree day sum (since freeze trigger) is lower than critical
                  ! value, then override the onset_flag set from soil water.

                  if (onset_gddflag(p) == 1._r8 .and. onset_gdd(p) < crit_onset_gdd) onset_flag(p) = 0._r8
               end if

               ! only allow onset if dayl > 6hrs
               if (onset_flag(p) == 1._r8 .and. dayl(g) <= secspqtrday) then
                  onset_flag(p) = 0._r8
               end if

               ! if this is the beginning of the onset period
               ! then reset the phenology flags and indices

               if (onset_flag(p) == 1._r8) then
                  dormant_flag(p) = 0._r8
                  days_active(p) = 0._r8
                  onset_gddflag(p) = 0._r8
                  onset_fdd(p) = 0._r8
                  onset_gdd(p) = 0._r8
                  onset_swi(p) = 0._r8
                  onset_counter(p) = ndays_on * secspday

                  ! call subroutine to move all the storage pools into transfer pools,
                  ! where they will be transfered to displayed growth over the onset period.
                  ! this code was originally handled with call cn_storage_to_xfer(p)
                  ! inlined during vectorization

                  ! set carbon fluxes for shifting storage pools to transfer pools
                  leafc_storage_to_xfer(p)  = fstor2tran * leafc_storage(p)/dt
                  frootc_storage_to_xfer(p) = fstor2tran * frootc_storage(p)/dt
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemc_storage_to_xfer(p)  = fstor2tran * livestemc_storage(p)/dt
                     deadstemc_storage_to_xfer(p)  = fstor2tran * deadstemc_storage(p)/dt
                     livecrootc_storage_to_xfer(p) = fstor2tran * livecrootc_storage(p)/dt
                     deadcrootc_storage_to_xfer(p) = fstor2tran * deadcrootc_storage(p)/dt
                     gresp_storage_to_xfer(p)      = fstor2tran * gresp_storage(p)/dt
                  end if

                  ! set nitrogen fluxes for shifting storage pools to transfer pools
                  leafn_storage_to_xfer(p)  = fstor2tran * leafn_storage(p)/dt
                  frootn_storage_to_xfer(p) = fstor2tran * frootn_storage(p)/dt
                  if (woody(ivt(p)) == 1.0_r8) then
                     livestemn_storage_to_xfer(p)  = fstor2tran * livestemn_storage(p)/dt
                     deadstemn_storage_to_xfer(p)  = fstor2tran * deadstemn_storage(p)/dt
                     livecrootn_storage_to_xfer(p) = fstor2tran * livecrootn_storage(p)/dt
                     deadcrootn_storage_to_xfer(p) = fstor2tran * deadcrootn_storage(p)/dt
                  end if
               end if

               ! test for switching from growth period to offset period
            else if (offset_flag(p) == 0._r8) then

               ! if soil water potential lower than critical value, accumulate
               ! as stress in offset soil water index

               if (psi <= soilpsi_off) then
                  offset_swi(p) = offset_swi(p) + fracday

                  ! if the offset soil water index exceeds critical value, and
                  ! if this is not the middle of a previously initiated onset period,
                  ! then set flag to start the offset period and reset index variables

                  if (offset_swi(p) >= crit_offset_swi .and. onset_flag(p) == 0._r8) offset_flag(p) = 1._r8

                  ! if soil water potential higher than critical value, reduce the
                  ! offset water stress index.  By this mechanism, there must be a
                  ! sustained period of water stress to initiate offset.

               else if (psi >= soilpsi_on) then
                  offset_swi(p) = offset_swi(p) - fracday
                  offset_swi(p) = max(offset_swi(p),0._r8)
               end if

               ! decrease freezing day accumulator for warm soil
               if (offset_fdd(p) > 0._r8 .and. soilt > SHR_CONST_TKFRZ) then
                  offset_fdd(p) = offset_fdd(p) - fracday
                  offset_fdd(p) = max(0._r8, offset_fdd(p))
               end if

               ! increase freezing day accumulator for cold soil
               if (soilt <= SHR_CONST_TKFRZ) then
                  offset_fdd(p) = offset_fdd(p) + fracday

                  ! if freezing degree day sum is greater than critical value, initiate offset
                  if (offset_fdd(p) > crit_offset_fdd .and. onset_flag(p) == 0._r8) offset_flag(p) = 1._r8
               end if

               ! force offset if daylength is < 6 hrs
               if (dayl(g) <= secspqtrday) then
                  offset_flag(p) = 1._r8
               end if

               ! if this is the beginning of the offset period
               ! then reset flags and indices
               if (offset_flag(p) == 1._r8) then
                  offset_fdd(p) = 0._r8
                  offset_swi(p) = 0._r8
                  offset_counter(p) = ndays_off * secspday
                  prev_leafc_to_litter(p) = 0._r8
                  prev_frootc_to_litter(p) = 0._r8
               end if
            end if

            ! keep track of number of days since last dormancy for control on
            ! fraction of new growth to send to storage for next growing season

            if (dormant_flag(p) == 0.0_r8) then
               days_active(p) = days_active(p) + fracday
            end if

            ! calculate long growing season factor (lgsf)
            ! only begin to calculate a lgsf greater than 0.0 once the number
            ! of days active exceeds days/year.
            lgsf(p) = max(min(3.0_r8*(days_active(p)-leaf_long(ivt(p))*dayspyr )/dayspyr, 1._r8),0._r8)
            ! RosieF. 5 Nov 2015.  Changed this such that the increase in leaf turnover is faster after
            ! trees enter the 'fake evergreen' state. Otherwise, they have a whole year of 
            ! cheating, with less litterfall than they should have, resulting in very high LAI. 
            ! Further, the 'fake evergreen' state (where lgsf>0) is entered at the end of a single leaf lifespan
            ! and not a whole year. The '3' is arbitrary, given that this entire system is quite abstract. 


            ! set background litterfall rate, when not in the phenological offset period
            if (offset_flag(p) == 1._r8) then
               bglfr(p) = 0._r8
            else
               ! calculate the background litterfall rate (bglfr)
               ! in units 1/s, based on leaf longevity (yrs) and correction for long growing season

               bglfr(p) = (1._r8/(leaf_long(ivt(p))*dayspyr*secspday))*lgsf(p)
            end if

            ! set background transfer rate when active but not in the phenological onset period
            if (onset_flag(p) == 1._r8) then
               bgtr(p) = 0._r8
            else
               ! the background transfer rate is calculated as the rate that would result
               ! in complete turnover of the storage pools in one year at steady state,
               ! once lgsf has reached 1.0 (after 730 days active).

               bgtr(p) = (1._r8/(dayspyr*secspday))*lgsf(p)

               ! set carbon fluxes for shifting storage pools to transfer pools

               ! reduced the amount of stored carbon flowing to display pool by only counting the delta
               ! between leafc and leafc_store in the flux. RosieF, Nov5 2015. 
               leafc_storage_to_xfer(p)  = max(0.0_r8,(leafc_storage(p)-leafc(p))) * bgtr(p)
               frootc_storage_to_xfer(p) = max(0.0_r8,(frootc_storage(p)-frootc(p))) * bgtr(p)
               if (woody(ivt(p)) == 1.0_r8) then
                  livestemc_storage_to_xfer(p)  = livestemc_storage(p) * bgtr(p)
                  deadstemc_storage_to_xfer(p)  = deadstemc_storage(p) * bgtr(p)
                  livecrootc_storage_to_xfer(p) = livecrootc_storage(p) * bgtr(p)
                  deadcrootc_storage_to_xfer(p) = deadcrootc_storage(p) * bgtr(p)
                  gresp_storage_to_xfer(p)      = gresp_storage(p) * bgtr(p)
               end if

               ! set nitrogen fluxes for shifting storage pools to transfer pools
               leafn_storage_to_xfer(p)  = leafn_storage(p) * bgtr(p)
               frootn_storage_to_xfer(p) = frootn_storage(p) * bgtr(p)
               if (woody(ivt(p)) == 1.0_r8) then
                  livestemn_storage_to_xfer(p)  = livestemn_storage(p) * bgtr(p)
                  deadstemn_storage_to_xfer(p)  = deadstemn_storage(p) * bgtr(p)
                  livecrootn_storage_to_xfer(p) = livecrootn_storage(p) * bgtr(p)
                  deadcrootn_storage_to_xfer(p) = deadcrootn_storage(p) * bgtr(p)
               end if
            end if

         end if ! end if stress deciduous

      end do ! end of patch loop

    end associate

  end subroutine CNStressDecidPhenology

  !-----------------------------------------------------------------------
  subroutine CropPhenology(num_pcropp, filter_pcropp                     , &
       waterdiagnosticbulk_inst, temperature_inst, crop_inst, canopystate_inst, cnveg_state_inst , &
       cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst,&
       c13_cnveg_carbonstate_inst, c14_cnveg_carbonstate_inst)

    ! !DESCRIPTION:
    ! Code from AgroIBIS to determine crop phenology and code from CN to
    ! handle CN fluxes during the phenological onset                       & offset periods.
    
    ! !USES:
    use clm_time_manager , only : get_curr_date, get_curr_calday, get_days_per_year, get_rad_step_size
    use pftconMod        , only : ntmp_corn, nswheat, nwwheat, ntmp_soybean
    use pftconMod        , only : nirrig_tmp_corn, nirrig_swheat, nirrig_wwheat, nirrig_tmp_soybean
    use pftconMod        , only : ntrp_corn, nsugarcane, ntrp_soybean, ncotton, nrice
    use pftconMod        , only : nirrig_trp_corn, nirrig_sugarcane, nirrig_trp_soybean
    use pftconMod        , only : nirrig_cotton, nirrig_rice
    use pftconMod        , only : nmiscanthus, nirrig_miscanthus, nswitchgrass, nirrig_switchgrass
    
    use clm_varcon       , only : spval, secspday
    use clm_varctl       , only : use_fertilizer 
    use clm_varctl       , only : use_c13, use_c14
    use clm_varcon       , only : c13ratio, c14ratio
    !
    ! !ARGUMENTS:
    integer                        , intent(in)    :: num_pcropp       ! number of prog crop patches in filter
    integer                        , intent(in)    :: filter_pcropp(:) ! filter for prognostic crop patches
    type(waterdiagnosticbulk_type)          , intent(in)    :: waterdiagnosticbulk_inst
    type(temperature_type)         , intent(in)    :: temperature_inst
    type(crop_type)                , intent(inout) :: crop_inst
    type(canopystate_type)         , intent(in)    :: canopystate_inst
    type(cnveg_state_type)         , intent(inout) :: cnveg_state_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type) , intent(inout) :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: c13_cnveg_carbonstate_inst
    type(cnveg_carbonstate_type)   , intent(inout) :: c14_cnveg_carbonstate_inst
    !
    ! LOCAL VARAIBLES:
    integer kyr       ! current year
    integer kmo       ! month of year  (1, ..., 12)
    integer kda       ! day of month   (1, ..., 31)
    integer mcsec     ! seconds of day (0, ..., seconds/day)
    integer jday      ! julian day of the year
    integer fp,p      ! patch indices
    integer c         ! column indices
    integer g         ! gridcell indices
    integer h         ! hemisphere indices
    integer idpp      ! number of days past planting
    real(r8) :: dtrad ! radiation time step delta t (seconds)
    real(r8) dayspyr  ! days per year
    real(r8) crmcorn  ! comparitive relative maturity for corn
    real(r8) ndays_on ! number of days to fertilize
    !------------------------------------------------------------------------

    associate(                                                                   & 
         ivt               =>    patch%itype                                     , & ! Input:  [integer  (:) ]  patch vegetation type                                
         
         leaf_long         =>    pftcon%leaf_long                              , & ! Input:  leaf longevity (yrs)                              
         leafcn            =>    pftcon%leafcn                                 , & ! Input:  leaf C:N (gC/gN)                                  
         manunitro         =>    pftcon%manunitro                              , & ! Input:  max manure to be applied in total (kgN/m2)
         mxmat             =>    pftcon%mxmat                                  , & ! Input:  
         minplanttemp      =>    pftcon%minplanttemp                           , & ! Input:  
         planttemp         =>    pftcon%planttemp                              , & ! Input:  
         gddmin            =>    pftcon%gddmin                                 , & ! Input:  
         hybgdd            =>    pftcon%hybgdd                                 , & ! Input:  
         lfemerg           =>    pftcon%lfemerg                                , & ! Input:  
         grnfill           =>    pftcon%grnfill                               , & ! Input:  

         t_ref2m_min       =>    temperature_inst%t_ref2m_min_patch            , & ! Input:  [real(r8) (:) ]  daily minimum of average 2 m height surface air temperature (K)
         t10               =>    temperature_inst%t_a10_patch                  , & ! Input:  [real(r8) (:) ]  10-day running mean of the 2 m temperature (K)    
         a5tmin            =>    temperature_inst%t_a5min_patch                , & ! Input:  [real(r8) (:) ]  5-day running mean of min 2-m temperature         
         a10tmin           =>    temperature_inst%t_a10min_patch               , & ! Input:  [real(r8) (:) ]  10-day running mean of min 2-m temperature        
         gdd020            =>    temperature_inst%gdd020_patch                 , & ! Input:  [real(r8) (:) ]  20 yr mean of gdd0                                
         gdd820            =>    temperature_inst%gdd820_patch                 , & ! Input:  [real(r8) (:) ]  20 yr mean of gdd8                                
         gdd1020           =>    temperature_inst%gdd1020_patch                , & ! Input:  [real(r8) (:) ]  20 yr mean of gdd10                               

         fertnitro         =>    crop_inst%fertnitro_patch                     , & ! Input:  [real(r8) (:) ]  fertilizer nitrogen
         hui               =>    crop_inst%gddplant_patch                      , & ! Input:  [real(r8) (:) ]  gdd since planting (gddplant)                    
         leafout           =>    crop_inst%gddtsoi_patch                       , & ! Input:  [real(r8) (:) ]  gdd from top soil layer temperature              
         harvdate          =>    crop_inst%harvdate_patch                      , & ! Output: [integer  (:) ]  harvest date                                       
         croplive          =>    crop_inst%croplive_patch                      , & ! Output: [logical  (:) ]  Flag, true if planted, not harvested               
         cropplant         =>    crop_inst%cropplant_patch                     , & ! Output: [logical  (:) ]  Flag, true if crop may be planted                  
         vf                =>    crop_inst%vf_patch                            , & ! Output: [real(r8) (:) ]  vernalization factor                              
         peaklai           =>  cnveg_state_inst%peaklai_patch                  , & ! Output: [integer  (:) ] 1: max allowed lai; 0: not at max                  
         tlai              =>    canopystate_inst%tlai_patch                   , & ! Input:  [real(r8) (:) ]  one-sided leaf area index, no burying by snow     
         
         idop              =>    cnveg_state_inst%idop_patch                   , & ! Output: [integer  (:) ]  date of planting                                   
         gddmaturity       =>    cnveg_state_inst%gddmaturity_patch            , & ! Output: [real(r8) (:) ]  gdd needed to harvest                             
         huileaf           =>    cnveg_state_inst%huileaf_patch                , & ! Output: [real(r8) (:) ]  heat unit index needed from planting to leaf emergence
         huigrain          =>    cnveg_state_inst%huigrain_patch               , & ! Output: [real(r8) (:) ]  same to reach vegetative maturity                 
         cumvd             =>    cnveg_state_inst%cumvd_patch                  , & ! Output: [real(r8) (:) ]  cumulative vernalization d?ependence?             
         hdidx             =>    cnveg_state_inst%hdidx_patch                  , & ! Output: [real(r8) (:) ]  cold hardening index?                             
         bglfr             =>    cnveg_state_inst%bglfr_patch                  , & ! Output: [real(r8) (:) ]  background litterfall rate (1/s)                  
         bgtr              =>    cnveg_state_inst%bgtr_patch                   , & ! Output: [real(r8) (:) ]  background transfer growth rate (1/s)             
         lgsf              =>    cnveg_state_inst%lgsf_patch                   , & ! Output: [real(r8) (:) ]  long growing season factor [0-1]                  
         onset_flag        =>    cnveg_state_inst%onset_flag_patch             , & ! Output: [real(r8) (:) ]  onset flag                                        
         offset_flag       =>    cnveg_state_inst%offset_flag_patch            , & ! Output: [real(r8) (:) ]  offset flag                                       
         onset_counter     =>    cnveg_state_inst%onset_counter_patch          , & ! Output: [real(r8) (:) ]  onset counter                                     
         offset_counter    =>    cnveg_state_inst%offset_counter_patch         , & ! Output: [real(r8) (:) ]  offset counter                                    
         
         leafc_xfer        =>    cnveg_carbonstate_inst%leafc_xfer_patch       , & ! Output: [real(r8) (:) ]  (gC/m2)   leaf C transfer                           

         crop_seedc_to_leaf =>   cnveg_carbonflux_inst%crop_seedc_to_leaf_patch, & ! Output: [real(r8) (:) ]  (gC/m2/s) seed source to leaf

         fert_counter      =>    cnveg_nitrogenflux_inst%fert_counter_patch    , & ! Output: [real(r8) (:) ]  >0 fertilize; <=0 not (seconds)                   
         leafn_xfer        =>    cnveg_nitrogenstate_inst%leafn_xfer_patch     , & ! Output: [real(r8) (:) ]  (gN/m2)   leaf N transfer                           
         crop_seedn_to_leaf =>   cnveg_nitrogenflux_inst%crop_seedn_to_leaf_patch, & ! Output: [real(r8) (:) ]  (gN/m2/s) seed source to leaf
         cphase            =>    crop_inst%cphase_patch                        , & ! Output: [real(r8) (:)]   phenology phase
         fert              =>    cnveg_nitrogenflux_inst%fert_patch              & ! Output: [real(r8) (:) ]  (gN/m2/s) fertilizer applied each timestep 
         )

      ! get time info
      dayspyr = get_days_per_year()
      jday    = get_curr_calday()
      call get_curr_date(kyr, kmo, kda, mcsec)
      dtrad   = real( get_rad_step_size(), r8 )

      if (use_fertilizer) then
       ndays_on = 20._r8 ! number of days to fertilize
      else
       ndays_on = 0._r8 ! number of days to fertilize
      end if

      do fp = 1, num_pcropp
         p = filter_pcropp(fp)
         c = patch%column(p)
         g = patch%gridcell(p)
         h = inhemi(p)

         ! background litterfall and transfer rates; long growing season factor

         bglfr(p) = 0._r8 ! this value changes later in a crop's life cycle
         bgtr(p)  = 0._r8
         lgsf(p)  = 0._r8

         ! ---------------------------------
         ! from AgroIBIS subroutine planting
         ! ---------------------------------

         ! in order to allow a crop to be planted only once each year
         ! initialize cropplant = .false., but hold it = .true. through the end of the year

         ! initialize other variables that are calculated for crops
         ! on an annual basis in cropresidue subroutine

         if ( jday == jdayyrstart(h) .and. mcsec == 0 )then

            ! make sure variables aren't changed at beginning of the year
            ! for a crop that is currently planted, such as
            ! WINTER TEMPERATE CEREAL = winter (wheat + barley + rye)
            ! represented here by the winter wheat pft

            if (.not. croplive(p))  then
               cropplant(p) = .false.
               idop(p)      = NOT_Planted

               ! keep next for continuous, annual winter temperate cereal crop;
               ! if we removed elseif,
               ! winter cereal grown continuously would amount to a cereal/fallow
               ! rotation because cereal would only be planted every other year

            else if (croplive(p) .and. (ivt(p) == nwwheat .or. ivt(p) == nirrig_wwheat)) then
               cropplant(p) = .false.
               !           else ! not possible to have croplive and ivt==cornORsoy? (slevis)
            end if

         end if

         if ( (.not. croplive(p)) .and. (.not. cropplant(p)) ) then

            ! gdd needed for * chosen crop and a likely hybrid (for that region) *
            ! to reach full physiological maturity

            ! based on accumulated seasonal average growing degree days from
            ! April 1 - Sept 30 (inclusive)
            ! for corn and soybeans in the United States -
            ! decided upon by what the typical average growing season length is
            ! and the gdd needed to reach maturity in those regions

            ! first choice is used for spring temperate cereal and/or soybeans and maize

            ! slevis: ibis reads xinpdate in io.f from control.crops.nc variable name 'plantdate'
            !         According to Chris Kucharik, the dataset of
            !         xinpdate was generated from a previous model run at 0.5 deg resolution

            ! winter temperate cereal : use gdd0 as a limit to plant winter cereal

            if (ivt(p) == nwwheat .or. ivt(p) == nirrig_wwheat) then

               ! add check to only plant winter cereal after other crops (soybean, maize)
               ! have been harvested

               ! *** remember order of planting is crucial - in terms of which crops you want
               ! to be grown in what order ***

               ! in this case, corn or soybeans are assumed to be planted before
               ! cereal would be in any particular year that both patches are allowed
               ! to grow in the same grid cell (e.g., double-cropping)

               ! slevis: harvdate below needs cropplant(p) above to be cropplant(p,ivt(p))
               !         where ivt(p) has rotated to winter cereal because
               !         cropplant through the end of the year for a harvested crop.
               !         Also harvdate(p) should be harvdate(p,ivt(p)) and should be
               !         updated on Jan 1st instead of at harvest (slevis)
               if (a5tmin(p)             /= spval                  .and. &
                    a5tmin(p)             <= minplanttemp(ivt(p))   .and. &
                    jday                  >= minplantjday(ivt(p),h) .and. &
                    (gdd020(p)            /= spval                  .and. &
                    gdd020(p)             >= gddmin(ivt(p)))) then

                  cumvd(p)       = 0._r8
                  hdidx(p)       = 0._r8
                  vf(p)          = 0._r8
                  croplive(p)    = .true.
                  cropplant(p)   = .true.
                  idop(p)        = jday
                  harvdate(p)    = NOT_Harvested
                  gddmaturity(p) = hybgdd(ivt(p))
                  leafc_xfer(p)  = initial_seed_at_planting
                  leafn_xfer(p)  = leafc_xfer(p) / leafcn(ivt(p)) ! with onset
                  crop_seedc_to_leaf(p) = leafc_xfer(p)/dt
                  crop_seedn_to_leaf(p) = leafn_xfer(p)/dt

                  ! because leafc_xfer is set above rather than incremneted through the normal process, must also set its isotope
                  ! pools here.  use totvegc_patch as the closest analogue if nonzero, and use initial value otherwise
                  if (use_c13) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c13_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c13ratio
                     endif
                  endif
                  if (use_c14) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c14_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c14ratio
                     endif
                  endif

                  ! latest possible date to plant winter cereal and after all other 
                  ! crops were harvested for that year

               else if (jday       >=  maxplantjday(ivt(p),h) .and. &
                    gdd020(p)  /= spval                   .and. &
                    gdd020(p)  >= gddmin(ivt(p))) then

                  cumvd(p)       = 0._r8
                  hdidx(p)       = 0._r8
                  vf(p)          = 0._r8
                  croplive(p)    = .true.
                  cropplant(p)   = .true.
                  idop(p)        = jday
                  harvdate(p)    = NOT_Harvested
                  gddmaturity(p) = hybgdd(ivt(p))
                  leafc_xfer(p)  = initial_seed_at_planting
                  leafn_xfer(p)  = leafc_xfer(p) / leafcn(ivt(p)) ! with onset
                  crop_seedc_to_leaf(p) = leafc_xfer(p)/dt
                  crop_seedn_to_leaf(p) = leafn_xfer(p)/dt

                  ! because leafc_xfer is set above rather than incremneted through the normal process, must also set its isotope
                  ! pools here.  use totvegc_patch as the closest analogue if nonzero, and use initial value otherwise
                  if (use_c13) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c13_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c13ratio
                     endif
                  endif
                  if (use_c14) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c14_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c14ratio
                     endif
                  endif
               else
                  gddmaturity(p) = 0._r8
               end if

            else ! not winter cereal... slevis: added distinction between NH and SH
               ! slevis: The idea is that jday will equal idop sooner or later in the year
               !         while the gdd part is either true or false for the year.
               if (t10(p) /= spval.and. a10tmin(p) /= spval   .and. &
                    t10(p)     > planttemp(ivt(p))             .and. &
                    a10tmin(p) > minplanttemp(ivt(p))          .and. &
                    jday       >= minplantjday(ivt(p),h)       .and. &
                    jday       <= maxplantjday(ivt(p),h)       .and. &
                    t10(p) /= spval .and. a10tmin(p) /= spval  .and. &
                    gdd820(p) /= spval                         .and. &
                    gdd820(p) >= gddmin(ivt(p))) then

                  ! impose limit on growing season length needed
                  ! for crop maturity - for cold weather constraints
                  croplive(p)  = .true.
                  cropplant(p) = .true.
                  idop(p)      = jday
                  harvdate(p)  = NOT_Harvested

                  ! go a specified amount of time before/after
                  ! climatological date
                  if (ivt(p) == ntmp_soybean .or. ivt(p) == nirrig_tmp_soybean .or. &
                       ivt(p) == ntrp_soybean .or. ivt(p) == nirrig_trp_soybean) then
                     gddmaturity(p) = min(gdd1020(p), hybgdd(ivt(p)))
                  end if
                  
                  if (ivt(p) == ntmp_corn .or. ivt(p) == nirrig_tmp_corn .or. &
                      ivt(p) == ntrp_corn .or. ivt(p) == nirrig_trp_corn .or. &
                      ivt(p) == nsugarcane .or. ivt(p) == nirrig_sugarcane .or. &
                      ivt(p) == nmiscanthus .or. ivt(p) == nirrig_miscanthus .or. &
                      ivt(p) == nswitchgrass .or. ivt(p) == nirrig_switchgrass) then
                     gddmaturity(p) = max(950._r8, min(gdd820(p)*0.85_r8, hybgdd(ivt(p))))
                     gddmaturity(p) = max(950._r8, min(gddmaturity(p)+150._r8, 1850._r8))
                  end if
                  if (ivt(p) == nswheat .or. ivt(p) == nirrig_swheat .or. &
                      ivt(p) == ncotton .or. ivt(p) == nirrig_cotton .or. &
                      ivt(p) == nrice   .or. ivt(p) == nirrig_rice) then
                     gddmaturity(p) = min(gdd020(p), hybgdd(ivt(p)))
                  end if

                  leafc_xfer(p)  = initial_seed_at_planting
                  leafn_xfer(p) = leafc_xfer(p) / leafcn(ivt(p)) ! with onset
                  crop_seedc_to_leaf(p) = leafc_xfer(p)/dt
                  crop_seedn_to_leaf(p) = leafn_xfer(p)/dt

                  ! because leafc_xfer is set above rather than incremneted through the normal process, must also set its isotope
                  ! pools here.  use totvegc_patch as the closest analogue if nonzero, and use initial value otherwise
                  if (use_c13) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c13_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c13ratio
                     endif
                  endif
                  if (use_c14) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c14_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c14ratio
                     endif
                  endif


                  ! If hit the max planting julian day -- go ahead and plant
               else if (jday == maxplantjday(ivt(p),h) .and. gdd820(p) > 0._r8 .and. &
                    gdd820(p) /= spval ) then
                  croplive(p)  = .true.
                  cropplant(p) = .true.
                  idop(p)      = jday
                  harvdate(p)  = NOT_Harvested

                  if (ivt(p) == ntmp_soybean .or. ivt(p) == nirrig_tmp_soybean .or. &
                      ivt(p) == ntrp_soybean .or. ivt(p) == nirrig_trp_soybean) then
                     gddmaturity(p) = min(gdd1020(p), hybgdd(ivt(p)))
                  end if
                  
                  if (ivt(p) == ntmp_corn .or. ivt(p) == nirrig_tmp_corn .or. &
                      ivt(p) == ntrp_corn .or. ivt(p) == nirrig_trp_corn .or. &
                      ivt(p) == nsugarcane .or. ivt(p) == nirrig_sugarcane .or. &
                      ivt(p) == nmiscanthus .or. ivt(p) == nirrig_miscanthus .or. &
                      ivt(p) == nswitchgrass .or. ivt(p) == nirrig_switchgrass) then
                     gddmaturity(p) = max(950._r8, min(gdd820(p)*0.85_r8, hybgdd(ivt(p))))
                  end if
                  if (ivt(p) == nswheat .or. ivt(p) == nirrig_swheat .or. &
                      ivt(p) == ncotton .or. ivt(p) == nirrig_cotton .or. &
                      ivt(p) == nrice   .or. ivt(p) == nirrig_rice) then
                     gddmaturity(p) = min(gdd020(p), hybgdd(ivt(p)))
                  end if

                  leafc_xfer(p)  = initial_seed_at_planting
                  leafn_xfer(p) = leafc_xfer(p) / leafcn(ivt(p)) ! with onset
                  crop_seedc_to_leaf(p) = leafc_xfer(p)/dt
                  crop_seedn_to_leaf(p) = leafn_xfer(p)/dt

                  ! because leafc_xfer is set above rather than incremneted through the normal process, must also set its isotope
                  ! pools here.  use totvegc_patch as the closest analogue if nonzero, and use initial value otherwise
                  if (use_c13) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c13_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c13ratio
                     endif
                  endif
                  if (use_c14) then
                     if ( cnveg_carbonstate_inst%totvegc_patch(p) .gt. 0._r8) then
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * &
                             c14_cnveg_carbonstate_inst%totvegc_patch(p) / cnveg_carbonstate_inst%totvegc_patch(p)
                     else
                        c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = leafc_xfer(p) * c14ratio
                     endif
                  endif

               else
                  gddmaturity(p) = 0._r8
               end if
            end if ! crop patch distinction

            ! crop phenology (gdd thresholds) controlled by gdd needed for
            ! maturity (physiological) which is based on the average gdd
            ! accumulation and hybrids in United States from April 1 - Sept 30

            ! calculate threshold from phase 1 to phase 2:
            ! threshold for attaining leaf emergence (based on fraction of
            ! gdd(i) -- climatological average)
            ! Hayhoe and Dwyer, 1990, Can. J. Soil Sci 70:493-497
            ! Carlson and Gage, 1989, Agric. For. Met., 45: 313-324
            ! J.T. Ritchie, 1991: Modeling Plant and Soil systems

            huileaf(p) = lfemerg(ivt(p)) * gddmaturity(p) ! 3-7% in cereal

            ! calculate threshhold from phase 2 to phase 3:
            ! from leaf emergence to beginning of grain-fill period
            ! this hypothetically occurs at the end of tassling, not the beginning
            ! tassel initiation typically begins at 0.5-0.55 * gddmaturity

            ! calculate linear relationship between huigrain fraction and relative
            ! maturity rating for maize

            if (ivt(p) == ntmp_corn .or. ivt(p) == nirrig_tmp_corn .or. &
                ivt(p) == ntrp_corn .or. ivt(p) == nirrig_trp_corn .or. &
                ivt(p) == nsugarcane .or. ivt(p) == nirrig_sugarcane .or. &
                ivt(p) == nmiscanthus .or. ivt(p) == nirrig_miscanthus .or. &
                ivt(p) == nswitchgrass .or. ivt(p) == nirrig_switchgrass) then
               ! the following estimation of crmcorn from gddmaturity is based on a linear
               ! regression using data from Pioneer-brand corn hybrids (Kucharik, 2003,
               ! Earth Interactions 7:1-33: fig. 2)
               crmcorn = max(73._r8, min(135._r8, (gddmaturity(p)+ 53.683_r8)/13.882_r8))

               ! the following adjustment of grnfill based on crmcorn is based on a tuning
               ! of Agro-IBIS to give reasonable results for max LAI and the seasonal
               ! progression of LAI growth (pers. comm. C. Kucharik June 10, 2010)
               huigrain(p) = -0.002_r8  * (crmcorn - 73._r8) + grnfill(ivt(p))

               huigrain(p) = min(max(huigrain(p), grnfill(ivt(p))-0.1_r8), grnfill(ivt(p)))
               huigrain(p) = huigrain(p) * gddmaturity(p)     ! Cabelguenne et
            else
               huigrain(p) = grnfill(ivt(p)) * gddmaturity(p) ! al. 1999
            end if

         end if ! crop not live nor planted

         ! ----------------------------------
         ! from AgroIBIS subroutine phenocrop
         ! ----------------------------------

         ! all of the phenology changes are based on the total number of gdd needed
         ! to change to the next phase - based on fractions of the total gdd typical
         ! for  that region based on the April 1 - Sept 30 window of development

         ! crop phenology (gdd thresholds) controlled by gdd needed for
         ! maturity (physiological) which is based on the average gdd
         ! accumulation and hybrids in United States from April 1 - Sept 30

         ! Phase 1: Planting to leaf emergence (now in CNAllocation)
         ! Phase 2: Leaf emergence to beginning of grain fill (general LAI accumulation)
         ! Phase 3: Grain fill to physiological maturity and harvest (LAI decline)
         ! Harvest: if gdd past grain fill initiation exceeds limit
         ! or number of days past planting reaches a maximum, the crop has
         ! reached physiological maturity and plant is harvested;
         ! crop could be live or dead at this stage - these limits
         ! could lead to reaching physiological maturity or determining
         ! a harvest date for a crop killed by an early frost (see next comments)
         ! --- --- ---
         ! keeping comments without the code (slevis):
         ! if minimum temperature, t_ref2m_min <= freeze kill threshold, tkill
         ! for 3 consecutive days and lai is above a minimum,
         ! plant will be damaged/killed. This function is more for spring freeze events
         ! or for early fall freeze events

         ! spring temperate cereal is affected by this, winter cereal kill function
         ! is determined in crops.f - is a more elaborate function of
         ! cold hardening of the plant

         ! currently simulates too many grid cells killed by freezing temperatures

         ! removed on March 12 2002 - C. Kucharik
         ! until it can be a bit more refined, or used at a smaller scale.
         ! we really have no way of validating this routine
         ! too difficult to implement on 0.5 degree scale grid cells
         ! --- --- ---

         onset_flag(p)  = 0._r8 ! CN terminology to trigger certain
         offset_flag(p) = 0._r8 ! carbon and nitrogen transfers

         if (croplive(p)) then
            cphase(p) = 1._r8

            ! call vernalization if winter temperate cereal planted, living, and the
            ! vernalization factor is not 1;
            ! vf affects the calculation of gddtsoi & gddplant

            if (t_ref2m_min(p) < 1.e30_r8 .and. vf(p) /= 1._r8 .and. &
               (ivt(p) == nwwheat .or. ivt(p) == nirrig_wwheat)) then
               call vernalization(p, &
                    canopystate_inst, temperature_inst, waterdiagnosticbulk_inst, cnveg_state_inst, &
                    crop_inst)
            end if

            ! days past planting may determine harvest

            if (jday >= idop(p)) then
               idpp = jday - idop(p)
            else
               idpp = int(dayspyr) + jday - idop(p)
            end if

            ! onset_counter initialized to zero when .not. croplive
            ! offset_counter relevant only at time step of harvest

            onset_counter(p) = onset_counter(p) - dt

            ! enter phase 2 onset for one time step:
            ! transfer seed carbon to leaf emergence

            if (peaklai(p) >= 1) then
               hui(p) = max(hui(p),huigrain(p))
            endif

            if (leafout(p) >= huileaf(p) .and. hui(p) < huigrain(p) .and. idpp < mxmat(ivt(p))) then
               cphase(p) = 2._r8
               if (abs(onset_counter(p)) > 1.e-6_r8) then
                  onset_flag(p)    = 1._r8
                  onset_counter(p) = dt
                    fert_counter(p)  = ndays_on * secspday
                    if (ndays_on .gt. 0) then
                       fert(p) = (manunitro(ivt(p)) * 1000._r8 + fertnitro(p))/ fert_counter(p)
                    else
                       fert(p) = 0._r8
                    end if
               else
                  ! this ensures no re-entry to onset of phase2
                  ! b/c onset_counter(p) = onset_counter(p) - dt
                  ! at every time step

                  onset_counter(p) = dt
               end if

               ! enter harvest for one time step:
               ! - transfer live biomass to litter and to crop yield
               ! - send xsmrpool to the atmosphere
               ! if onset and harvest needed to last longer than one timestep
               ! the onset_counter would change from dt and you'd need to make
               ! changes to the offset subroutine below

            else if (hui(p) >= gddmaturity(p) .or. idpp >= mxmat(ivt(p))) then
               if (harvdate(p) >= NOT_Harvested) harvdate(p) = jday
               croplive(p) = .false.     ! no re-entry in greater if-block
               cphase(p) = 4._r8
               if (tlai(p) > 0._r8) then ! plant had emerged before harvest
                  offset_flag(p) = 1._r8
                  offset_counter(p) = dt
               else                      ! plant never emerged from the ground
                  ! Revert planting transfers; this will replenish the crop seed deficit.
                  ! We subtract from any existing value in crop_seedc_to_leaf /
                  ! crop_seedn_to_leaf in the unlikely event that we enter this block of
                  ! code in the same time step where the planting transfer originally
                  ! occurred.
                  crop_seedc_to_leaf(p) = crop_seedc_to_leaf(p) - leafc_xfer(p)/dt
                  crop_seedn_to_leaf(p) = crop_seedn_to_leaf(p) - leafn_xfer(p)/dt
                  leafc_xfer(p) = 0._r8
                  leafn_xfer(p) = leafc_xfer(p) / leafcn(ivt(p))
                  if (use_c13) then
                     c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = 0._r8
                  endif
                  if (use_c14) then
                     c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = 0._r8
                  endif

               end if

               ! enter phase 3 while previous criteria fail and next is true;
               ! in terms of order, phase 3 occurs before harvest, but when
               ! harvest *can* occur, we want it to have first priority.
               ! AgroIBIS uses a complex formula for lai decline.
               ! Use CN's simple formula at least as a place holder (slevis)

            else if (hui(p) >= huigrain(p)) then
               cphase(p) = 3._r8
               bglfr(p) = 1._r8/(leaf_long(ivt(p))*dayspyr*secspday)
            end if

            ! continue fertilizer application while in phase 2;
            ! assumes that onset of phase 2 took one time step only

              if (fert_counter(p) <= 0._r8) then
                 fert(p) = 0._r8
              else ! continue same fert application every timestep
                 fert_counter(p) = fert_counter(p) - dtrad
              end if

         else   ! crop not live
            ! next 2 lines conserve mass if leaf*_xfer > 0 due to interpinic.
            ! We subtract from any existing value in crop_seedc_to_leaf /
            ! crop_seedn_to_leaf in the unlikely event that we enter this block of
            ! code in the same time step where the planting transfer originally
            ! occurred.
            crop_seedc_to_leaf(p) = crop_seedc_to_leaf(p) - leafc_xfer(p)/dt
            crop_seedn_to_leaf(p) = crop_seedn_to_leaf(p) - leafn_xfer(p)/dt
            onset_counter(p) = 0._r8
            leafc_xfer(p) = 0._r8
            leafn_xfer(p) = leafc_xfer(p) / leafcn(ivt(p))
            if (use_c13) then
               c13_cnveg_carbonstate_inst%leafc_xfer_patch(p) = 0._r8
            endif
            if (use_c14) then
               c14_cnveg_carbonstate_inst%leafc_xfer_patch(p) = 0._r8
            endif
         end if ! croplive

      end do ! prognostic crops loop

    end associate

  end subroutine CropPhenology

  !-----------------------------------------------------------------------
  subroutine CropPhenologyInit(bounds)
    !
    ! !DESCRIPTION:
    ! Initialization of CropPhenology. Must be called after time-manager is
    ! initialized, and after pftcon file is read in.
    !
    ! !USES:
    use pftconMod       , only: npcropmin, npcropmax
    use clm_time_manager, only: get_calday
    !
    ! !ARGUMENTS:
    type(bounds_type), intent(in) :: bounds  
    !
    ! LOCAL VARAIBLES:
    integer           :: p,g,n,i                     ! indices
    !------------------------------------------------------------------------

    allocate( inhemi(bounds%begp:bounds%endp) )

    allocate( minplantjday(0:maxveg,inSH)) ! minimum planting julian day
    allocate( maxplantjday(0:maxveg,inSH)) ! minimum planting julian day

    ! Julian day for the start of the year (mid-winter)
    jdayyrstart(inNH) =   1
    jdayyrstart(inSH) = 182

    ! Convert planting dates into julian day
    minplantjday(:,:) = huge(1)
    maxplantjday(:,:) = huge(1)
    do n = npcropmin, npcropmax
       if (pftcon%is_pft_known_to_model(n)) then
          minplantjday(n, inNH) = int( get_calday( pftcon%mnNHplantdate(n), 0 ) )
          maxplantjday(n, inNH) = int( get_calday( pftcon%mxNHplantdate(n), 0 ) )

          minplantjday(n, inSH) = int( get_calday( pftcon%mnSHplantdate(n), 0 ) )
          maxplantjday(n, inSH) = int( get_calday( pftcon%mxSHplantdate(n), 0 ) )
       end if
    end do

    ! Figure out what hemisphere each PATCH is in
    do p = bounds%begp, bounds%endp
       g = patch%gridcell(p)
       ! Northern hemisphere
       if ( grc%latdeg(g) > 0.0_r8 )then
          inhemi(p) = inNH
       else
          inhemi(p) = inSH
       end if
    end do

    !
    ! Constants for Crop vernalization
    !
    ! photoperiod factor calculation
    ! genetic constant - can be modified

    p1d = 0.004_r8  ! average for genotypes from Ritchey, 1991.
    ! Modeling plant & soil systems: Wheat phasic developmt
    p1v = 0.003_r8  ! average for genotypes from Ritchey, 1991.

    hti   = 1._r8
    tbase = 0._r8

  end subroutine CropPhenologyInit

  !-----------------------------------------------------------------------
  subroutine vernalization(p, &
       canopystate_inst, temperature_inst, waterdiagnosticbulk_inst, cnveg_state_inst, crop_inst)
    !
    ! !DESCRIPTION:
    !
    ! * * * only call for winter temperate cereal * * *
    !
    ! subroutine calculates vernalization and photoperiod effects on
    ! gdd accumulation in winter temperate cereal varieties. Thermal time accumulation
    ! is reduced in 1st period until plant is fully vernalized. During this
    ! time of emergence to spikelet formation, photoperiod can also have a
    ! drastic effect on plant development.
    !
    ! !ARGUMENTS:
    integer                , intent(in)    :: p    ! PATCH index running over
    type(canopystate_type) , intent(in)    :: canopystate_inst
    type(temperature_type) , intent(in)    :: temperature_inst
    type(waterdiagnosticbulk_type)  , intent(in)    :: waterdiagnosticbulk_inst
    type(cnveg_state_type) , intent(inout) :: cnveg_state_inst
    type(crop_type)        , intent(inout) :: crop_inst
    !
    ! LOCAL VARAIBLES:
    real(r8) tcrown                     ! ?
    real(r8) vd, vd1, vd2               ! vernalization dependence
    real(r8) tkil                       ! Freeze kill threshold
    integer  c,g                        ! indices
    !------------------------------------------------------------------------

    associate(                                               & 
         tlai        => canopystate_inst%tlai_patch        , & ! Input:  [real(r8) (:) ]  one-sided leaf area index, no burying by snow     

         t_ref2m     => temperature_inst%t_ref2m_patch     , & ! Input:  [real(r8) (:) ]  2 m height surface air temperature (K)            
         t_ref2m_min => temperature_inst%t_ref2m_min_patch , & ! Input:  [real(r8) (:) ] daily minimum of average 2 m height surface air temperature (K)
         t_ref2m_max => temperature_inst%t_ref2m_max_patch , & ! Input:  [real(r8) (:) ] daily maximum of average 2 m height surface air temperature (K)

         snow_depth  => waterdiagnosticbulk_inst%snow_depth_col     , & ! Input:  [real(r8) (:) ]  snow height (m)                                   

         hdidx       => cnveg_state_inst%hdidx_patch       , & ! Output: [real(r8) (:) ]  cold hardening index?                             
         cumvd       => cnveg_state_inst%cumvd_patch       , & ! Output: [real(r8) (:) ]  cumulative vernalization d?ependence?             
         gddmaturity => cnveg_state_inst%gddmaturity_patch , & ! Output: [real(r8) (:) ]  gdd needed to harvest                             
         huigrain    => cnveg_state_inst%huigrain_patch    , & ! Output: [real(r8) (:) ]  heat unit index needed to reach vegetative maturity

         vf          => crop_inst%vf_patch                   & ! Output: [real(r8) (:) ]  vernalization factor for cereal
         )

      c = patch%column(p)

      ! for all equations - temperatures must be in degrees (C)
      ! calculate temperature of crown of crop (e.g., 3 cm soil temperature)
      ! snow depth in centimeters
      
      if (t_ref2m(p) < tfrz) then !slevis: t_ref2m inst of td=daily avg (K)
         tcrown = 2._r8 + (t_ref2m(p) - tfrz) * (0.4_r8 + 0.0018_r8 * &
              (min(snow_depth(c)*100._r8, 15._r8) - 15._r8)**2)
      else !slevis: snow_depth inst of adsnod=daily average (m)
         tcrown = t_ref2m(p) - tfrz
      end if

      ! vernalization factor calculation
      ! if vf(p) = 1.  then plant is fully vernalized - and thermal time
      ! accumulation in phase 1 will be unaffected
      ! refers to gddtsoi & gddplant, defined in the accumulation routines (slevis)
      ! reset vf, cumvd, and hdidx to 0 at planting of crop (slevis)

      if (t_ref2m_max(p) > tfrz) then
         if (t_ref2m_min(p) <= tfrz+15._r8) then
            vd1      = 1.4_r8 - 0.0778_r8 * tcrown
            vd2      = 0.5_r8 + 13.44_r8 / ((t_ref2m_max(p)-t_ref2m_min(p)+3._r8)**2) * tcrown
            vd       = max(0._r8, min(1._r8, vd1, vd2))
            cumvd(p) = cumvd(p) + vd
         end if

         if (cumvd(p) < 10._r8 .and. t_ref2m_max(p) > tfrz+30._r8) then
            cumvd(p) = cumvd(p) - 0.5_r8 * (t_ref2m_max(p) - tfrz - 30._r8)
         end if
         cumvd(p) = max(0._r8, cumvd(p))       ! must be > 0

         vf(p) = 1._r8 - p1v * (50._r8 - cumvd(p))
         vf(p) = max(0._r8, min(vf(p), 1._r8)) ! must be between 0 - 1
      end if

      ! calculate cold hardening of plant
      ! determines for winter cereal varieties whether the plant has completed
      ! a period of cold hardening to protect it from freezing temperatures. If
      ! not, then exposure could result in death or killing of plants.

      ! there are two distinct phases of hardening

      if (t_ref2m_min(p) <= tfrz-3._r8 .or. hdidx(p) /= 0._r8) then
         if (hdidx(p) >= hti) then   ! done with phase 1
            hdidx(p) = hdidx(p) + 0.083_r8
            hdidx(p) = min(hdidx(p), hti*2._r8)
         end if

         if (t_ref2m_max(p) >= tbase + tfrz + 10._r8) then
            hdidx(p) = hdidx(p) - 0.02_r8 * (t_ref2m_max(p)-tbase-tfrz-10._r8)
            if (hdidx(p) > hti) hdidx(p) = hdidx(p) - 0.02_r8 * (t_ref2m_max(p)-tbase-tfrz-10._r8)
            hdidx(p) = max(0._r8, hdidx(p))
         end if

      else if (tcrown >= tbase-1._r8) then
         if (tcrown <= tbase+8._r8) then
            hdidx(p) = hdidx(p) + 0.1_r8 - (tcrown-tbase+3.5_r8)**2 / 506._r8
            if (hdidx(p) >= hti .and. tcrown <= tbase + 0._r8) then
               hdidx(p) = hdidx(p) + 0.083_r8
               hdidx(p) = min(hdidx(p), hti*2._r8)
            end if
         end if

         if (t_ref2m_max(p) >= tbase + tfrz + 10._r8) then
            hdidx(p) = hdidx(p) - 0.02_r8 * (t_ref2m_max(p)-tbase-tfrz-10._r8)
            if (hdidx(p) > hti) hdidx(p) = hdidx(p) - 0.02_r8 * (t_ref2m_max(p)-tbase-tfrz-10._r8)
            hdidx(p) = max(0._r8, hdidx(p))
         end if
      end if

      ! calculate what the cereal killing temperature
      ! there is a linear inverse relationship between
      ! hardening of the plant and the killing temperature or
      ! threshold that the plant can withstand
      ! when plant is fully-hardened (hdidx = 2), the killing threshold is -18 C

      ! will have to develop some type of relationship that reduces LAI and
      ! biomass pools in response to cold damaged crop

      if (t_ref2m_min(p) <= tfrz - 6._r8) then
         tkil = (tbase - 6._r8) - 6._r8 * hdidx(p)
         if (tkil >= tcrown) then
            if ((0.95_r8 - 0.02_r8 * (tcrown - tkil)**2) >= 0.02_r8) then
               write (iulog,*)  'crop damaged by cold temperatures at p,c =', p,c
            else if (tlai(p) > 0._r8) then ! slevis: kill if past phase1
               gddmaturity(p) = 0._r8      !         by forcing through
               huigrain(p)    = 0._r8      !         harvest
               write (iulog,*)  '95% of crop killed by cold temperatures at p,c =', p,c
            end if
         end if
      end if

    end associate 

  end subroutine vernalization

  !-----------------------------------------------------------------------
  subroutine CNOnsetGrowth (num_soilp, filter_soilp, &
       cnveg_state_inst, &
       cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! Determines the flux of stored C and N from transfer pools to display
    ! pools during the phenological onset period.
    !
    ! !ARGUMENTS:
    integer                        , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                        , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(cnveg_state_type)             , intent(in)    :: cnveg_state_inst
    type(cnveg_carbonstate_type)   , intent(in)    :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type) , intent(in)    :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    integer :: p            ! indices
    integer :: fp           ! lake filter patch index
    real(r8):: t1           ! temporary variable
    !-----------------------------------------------------------------------

    associate(                                                                                             & 
         ivt                                 =>    patch%itype                                                   , & ! Input:  [integer   (:) ]  patch vegetation type                                

         woody                               =>    pftcon%woody                                                , & ! Input:  binary flag for woody lifeform (1=woody, 0=not woody)
         
         onset_flag                          =>    cnveg_state_inst%onset_flag_patch                           , & ! Input:  [real(r8)  (:) ]  onset flag                                        
         onset_counter                       =>    cnveg_state_inst%onset_counter_patch                        , & ! Input:  [real(r8)  (:) ]  onset days counter                                
         bgtr                                =>    cnveg_state_inst%bgtr_patch                                 , & ! Input:  [real(r8)  (:) ]  background transfer growth rate (1/s)             
         
         leafc_xfer                          =>    cnveg_carbonstate_inst%leafc_xfer_patch                     , & ! Input:  [real(r8)  (:) ]  (gC/m2) leaf C transfer                           
         frootc_xfer                         =>    cnveg_carbonstate_inst%frootc_xfer_patch                    , & ! Input:  [real(r8)  (:) ]  (gC/m2) fine root C transfer                      
         livestemc_xfer                      =>    cnveg_carbonstate_inst%livestemc_xfer_patch                 , & ! Input:  [real(r8)  (:) ]  (gC/m2) live stem C transfer                      
         deadstemc_xfer                      =>    cnveg_carbonstate_inst%deadstemc_xfer_patch                 , & ! Input:  [real(r8)  (:) ]  (gC/m2) dead stem C transfer                      
         livecrootc_xfer                     =>    cnveg_carbonstate_inst%livecrootc_xfer_patch                , & ! Input:  [real(r8)  (:) ]  (gC/m2) live coarse root C transfer               
         deadcrootc_xfer                     =>    cnveg_carbonstate_inst%deadcrootc_xfer_patch                , & ! Input:  [real(r8)  (:) ]  (gC/m2) dead coarse root C transfer               
         
         leafn_xfer                          =>    cnveg_nitrogenstate_inst%leafn_xfer_patch                   , & ! Input:  [real(r8)  (:) ]  (gN/m2) leaf N transfer                           
         frootn_xfer                         =>    cnveg_nitrogenstate_inst%frootn_xfer_patch                  , & ! Input:  [real(r8)  (:) ]  (gN/m2) fine root N transfer                      
         livestemn_xfer                      =>    cnveg_nitrogenstate_inst%livestemn_xfer_patch               , & ! Input:  [real(r8)  (:) ]  (gN/m2) live stem N transfer                      
         deadstemn_xfer                      =>    cnveg_nitrogenstate_inst%deadstemn_xfer_patch               , & ! Input:  [real(r8)  (:) ]  (gN/m2) dead stem N transfer                      
         livecrootn_xfer                     =>    cnveg_nitrogenstate_inst%livecrootn_xfer_patch              , & ! Input:  [real(r8)  (:) ]  (gN/m2) live coarse root N transfer               
         deadcrootn_xfer                     =>    cnveg_nitrogenstate_inst%deadcrootn_xfer_patch              , & ! Input:  [real(r8)  (:) ]  (gN/m2) dead coarse root N transfer               
         
         leafc_xfer_to_leafc                 =>    cnveg_carbonflux_inst%leafc_xfer_to_leafc_patch             , & ! Output:  [real(r8) (:) ]                                                    
         frootc_xfer_to_frootc               =>    cnveg_carbonflux_inst%frootc_xfer_to_frootc_patch           , & ! Output:  [real(r8) (:) ]                                                    
         livestemc_xfer_to_livestemc         =>    cnveg_carbonflux_inst%livestemc_xfer_to_livestemc_patch     , & ! Output:  [real(r8) (:) ]                                                    
         deadstemc_xfer_to_deadstemc         =>    cnveg_carbonflux_inst%deadstemc_xfer_to_deadstemc_patch     , & ! Output:  [real(r8) (:) ]                                                    
         livecrootc_xfer_to_livecrootc       =>    cnveg_carbonflux_inst%livecrootc_xfer_to_livecrootc_patch   , & ! Output:  [real(r8) (:) ]                                                    
         deadcrootc_xfer_to_deadcrootc       =>    cnveg_carbonflux_inst%deadcrootc_xfer_to_deadcrootc_patch   , & ! Output:  [real(r8) (:) ]                                                    
         
         leafn_xfer_to_leafn                 =>    cnveg_nitrogenflux_inst%leafn_xfer_to_leafn_patch           , & ! Output:  [real(r8) (:) ]                                                    
         frootn_xfer_to_frootn               =>    cnveg_nitrogenflux_inst%frootn_xfer_to_frootn_patch         , & ! Output:  [real(r8) (:) ]                                                    
         livestemn_xfer_to_livestemn         =>    cnveg_nitrogenflux_inst%livestemn_xfer_to_livestemn_patch   , & ! Output:  [real(r8) (:) ]                                                    
         deadstemn_xfer_to_deadstemn         =>    cnveg_nitrogenflux_inst%deadstemn_xfer_to_deadstemn_patch   , & ! Output:  [real(r8) (:) ]                                                    
         livecrootn_xfer_to_livecrootn       =>    cnveg_nitrogenflux_inst%livecrootn_xfer_to_livecrootn_patch , & ! Output:  [real(r8) (:) ]                                                    
         deadcrootn_xfer_to_deadcrootn       =>    cnveg_nitrogenflux_inst%deadcrootn_xfer_to_deadcrootn_patch   & ! Output:  [real(r8) (:) ]                                                    
         )

      ! patch loop
      do fp = 1,num_soilp
         p = filter_soilp(fp)

         ! only calculate these fluxes during onset period
         if (onset_flag(p) == 1._r8) then

            ! The transfer rate is a linearly decreasing function of time,
            ! going to zero on the last timestep of the onset period

            if (onset_counter(p) == dt) then
               t1 = 1.0_r8 / dt
            else
               t1 = 2.0_r8 / (onset_counter(p))
            end if
            leafc_xfer_to_leafc(p)   = t1 * leafc_xfer(p)
            frootc_xfer_to_frootc(p) = t1 * frootc_xfer(p)
            leafn_xfer_to_leafn(p)   = t1 * leafn_xfer(p)
            frootn_xfer_to_frootn(p) = t1 * frootn_xfer(p)
            if (woody(ivt(p)) == 1.0_r8) then
               livestemc_xfer_to_livestemc(p)   = t1 * livestemc_xfer(p)
               deadstemc_xfer_to_deadstemc(p)   = t1 * deadstemc_xfer(p)
               livecrootc_xfer_to_livecrootc(p) = t1 * livecrootc_xfer(p)
               deadcrootc_xfer_to_deadcrootc(p) = t1 * deadcrootc_xfer(p)
               livestemn_xfer_to_livestemn(p)   = t1 * livestemn_xfer(p)
               deadstemn_xfer_to_deadstemn(p)   = t1 * deadstemn_xfer(p)
               livecrootn_xfer_to_livecrootn(p) = t1 * livecrootn_xfer(p)
               deadcrootn_xfer_to_deadcrootn(p) = t1 * deadcrootn_xfer(p)
            end if

         end if ! end if onset period

         ! calculate the background rate of transfer growth (used for stress
         ! deciduous algorithm). In this case, all of the mass in the transfer
         ! pools should be moved to displayed growth in each timestep.

         if (bgtr(p) > 0._r8) then
            leafc_xfer_to_leafc(p)   = leafc_xfer(p) / dt
            frootc_xfer_to_frootc(p) = frootc_xfer(p) / dt
            leafn_xfer_to_leafn(p)   = leafn_xfer(p) / dt
            frootn_xfer_to_frootn(p) = frootn_xfer(p) / dt
            if (woody(ivt(p)) == 1.0_r8) then
               livestemc_xfer_to_livestemc(p)   = livestemc_xfer(p) / dt
               deadstemc_xfer_to_deadstemc(p)   = deadstemc_xfer(p) / dt
               livecrootc_xfer_to_livecrootc(p) = livecrootc_xfer(p) / dt
               deadcrootc_xfer_to_deadcrootc(p) = deadcrootc_xfer(p) / dt
               livestemn_xfer_to_livestemn(p)   = livestemn_xfer(p) / dt
               deadstemn_xfer_to_deadstemn(p)   = deadstemn_xfer(p) / dt
               livecrootn_xfer_to_livecrootn(p) = livecrootn_xfer(p) / dt
               deadcrootn_xfer_to_deadcrootn(p) = deadcrootn_xfer(p) / dt
            end if
         end if ! end if bgtr

      end do ! end patch loop

    end associate

  end subroutine CNOnsetGrowth

  !-----------------------------------------------------------------------
  subroutine CNOffsetLitterfall (num_soilp, filter_soilp, &
       cnveg_state_inst, cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! Determines the flux of C and N from displayed pools to litter
    ! pools during the phenological offset period.
    !
    ! !USES:
    use pftconMod        , only : npcropmin
    use pftconMod        , only : nmiscanthus, nirrig_miscanthus, nswitchgrass, nirrig_switchgrass
    
    use CNSharedParamsMod, only : use_fun
    use clm_varctl       , only : CNratio_floating    
    !
    ! !ARGUMENTS:
    integer                       , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                       , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(cnveg_state_type)        , intent(inout) :: cnveg_state_inst
    type(cnveg_carbonstate_type)  , intent(in)    :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type), intent(in)    :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)   , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type) , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    integer :: p, c         ! indices
    integer :: fp           ! lake filter patch index
    real(r8):: t1           ! temporary variable
    real(r8):: denom        ! temporary variable for divisor
    real(r8) :: ntovr_leaf  
    real(r8) :: fr_leafn_to_litter ! fraction of the nitrogen turnover that goes to litter; remaining fraction is retranslocated
    !-----------------------------------------------------------------------

    associate(                                                                           & 
         ivt                   =>    patch%itype                                       , & ! Input:  [integer  (:) ]  patch vegetation type                                

         leafcn                =>    pftcon%leafcn                                     , & ! Input:  leaf C:N (gC/gN) 
         
         biofuel_harvfrac      =>    pftcon%biofuel_harvfrac                           , & ! Input:  cut a fraction of leaf & stem for biofuel (-) 
                                          
         lflitcn               =>    pftcon%lflitcn                                    , & ! Input:  leaf litter C:N (gC/gN)                           
         frootcn               =>    pftcon%frootcn                                    , & ! Input:  fine root C:N (gC/gN)                             
         graincn               =>    pftcon%graincn                                    , & ! Input:  grain C:N (gC/gN)                                 

         offset_flag           =>    cnveg_state_inst%offset_flag_patch                , & ! Input:  [real(r8) (:) ]  offset flag                                       
         offset_counter        =>    cnveg_state_inst%offset_counter_patch             , & ! Input:  [real(r8) (:) ]  offset days counter                               

         leafc                 =>    cnveg_carbonstate_inst%leafc_patch                , & ! Input:  [real(r8) (:) ]  (gC/m2) leaf C                                    
         frootc                =>    cnveg_carbonstate_inst%frootc_patch               , & ! Input:  [real(r8) (:) ]  (gC/m2) fine root C                               
         grainc                =>    cnveg_carbonstate_inst%grainc_patch               , & ! Input:  [real(r8) (:) ]  (gC/m2) grain C                                   
         cropseedc_deficit     =>    cnveg_carbonstate_inst%cropseedc_deficit_patch    , & ! Input:  [real(r8) (:) ]  (gC/m2) crop seed C deficit
         livestemc             =>    cnveg_carbonstate_inst%livestemc_patch            , & ! Input:  [real(r8) (:) ]  (gC/m2) livestem C                                
         cropseedn_deficit     =>    cnveg_nitrogenstate_inst%cropseedn_deficit_patch  , & ! Input:  [real(r8) (:) ]  (gC/m2) crop seed N deficit
         livestemn             =>    cnveg_nitrogenstate_inst%livestemn_patch          , & ! Input:  [real(r8) (:) ]  (gN/m2) livestem N

         cpool_to_grainc       =>    cnveg_carbonflux_inst%cpool_to_grainc_patch       , & ! Input:  [real(r8) (:) ]  allocation to grain C (gC/m2/s)                   
         npool_to_grainn       =>    cnveg_nitrogenflux_inst%npool_to_grainn_patch     , & ! Input: [real(r8) (:)  ]  allocation to grain N (gN/m2/s)
         grainn                =>    cnveg_nitrogenstate_inst%grainn_patch             , & ! Input: [real(r8) (:)  ]  (kgN/m2) grain N
         cpool_to_livestemc    =>    cnveg_carbonflux_inst%cpool_to_livestemc_patch    , & ! Input:  [real(r8) (:) ]  allocation to live stem C (gC/m2/s)               
         cpool_to_leafc        =>    cnveg_carbonflux_inst%cpool_to_leafc_patch        , & ! Input:  [real(r8) (:) ]  allocation to leaf C (gC/m2/s)                    
         cpool_to_frootc       =>    cnveg_carbonflux_inst%cpool_to_frootc_patch       , & ! Input:  [real(r8) (:) ]  allocation to fine root C (gC/m2/s)               
         prev_leafc_to_litter  =>    cnveg_carbonflux_inst%prev_leafc_to_litter_patch  , & ! Output: [real(r8) (:) ]  previous timestep leaf C litterfall flux (gC/m2/s)
         prev_frootc_to_litter =>    cnveg_carbonflux_inst%prev_frootc_to_litter_patch , & ! Output: [real(r8) (:) ]  previous timestep froot C litterfall flux (gC/m2/s)
         leafc_to_litter       =>    cnveg_carbonflux_inst%leafc_to_litter_patch       , & ! Output: [real(r8) (:) ]  leaf C litterfall (gC/m2/s)                       
         frootc_to_litter      =>    cnveg_carbonflux_inst%frootc_to_litter_patch      , & ! Output: [real(r8) (:) ]  fine root C litterfall (gC/m2/s)                  
         livestemc_to_litter   =>    cnveg_carbonflux_inst%livestemc_to_litter_patch   , & ! Output: [real(r8) (:) ]  live stem C litterfall (gC/m2/s)                  
         grainc_to_food        =>    cnveg_carbonflux_inst%grainc_to_food_patch        , & ! Output: [real(r8) (:) ]  grain C to food (gC/m2/s)             
         grainc_to_seed        =>    cnveg_carbonflux_inst%grainc_to_seed_patch        , & ! Output: [real(r8) (:) ]  grain C to seed (gC/m2/s)
         leafc_to_biofuelc     =>    cnveg_carbonflux_inst%leafc_to_biofuelc_patch     , & ! Output: [real(r8) (:) ]  leaf C to biofuel C (gC/m2/s)
         livestemc_to_biofuelc =>    cnveg_carbonflux_inst%livestemc_to_biofuelc_patch , & ! Output: [real(r8) (:) ]  livestem C to biofuel C (gC/m2/s)
         leafn                 =>    cnveg_nitrogenstate_inst%leafn_patch              , & ! Input:  [real(r8) (:) ]  (gN/m2) leaf N      
         frootn                =>    cnveg_nitrogenstate_inst%frootn_patch             , & ! Input:  [real(r8) (:) ]  (gN/m2) fine root N                        

         livestemn_to_litter   =>    cnveg_nitrogenflux_inst%livestemn_to_litter_patch , & ! Output: [real(r8) (:) ]  livestem N to litter (gN/m2/s)                    
         grainn_to_food        =>    cnveg_nitrogenflux_inst%grainn_to_food_patch      , & ! Output: [real(r8) (:) ]  grain N to food (gN/m2/s)                                   
         grainn_to_seed        =>    cnveg_nitrogenflux_inst%grainn_to_seed_patch      , & ! Output: [real(r8) (:) ]  grain N to seed (gN/m2/s)
         leafn_to_biofueln     =>    cnveg_nitrogenflux_inst%leafn_to_biofueln_patch   , & ! Output: [real(r8) (:) ]  leaf N to biofuel N (gN/m2/s)
         livestemn_to_biofueln =>    cnveg_nitrogenflux_inst%livestemn_to_biofueln_patch, & ! Output: [real(r8) (:) ]  livestem N to biofuel N (gN/m2/s)     
         leafn_to_litter       =>    cnveg_nitrogenflux_inst%leafn_to_litter_patch     , & ! Output: [real(r8) (:) ]  leaf N litterfall (gN/m2/s)                       
         leafn_to_retransn     =>    cnveg_nitrogenflux_inst%leafn_to_retransn_patch   , & ! Input: [real(r8) (:) ]  leaf N to retranslocated N pool (gN/m2/s)         
         free_retransn_to_npool=>    cnveg_nitrogenflux_inst%free_retransn_to_npool_patch  , & ! Input: [real(r8) (:) ] free leaf N to retranslocated N pool (gN/m2/s)          
         paid_retransn_to_npool=>    cnveg_nitrogenflux_inst%retransn_to_npool_patch, & ! Input: [real(r8) (:) ] free leaf N to retranslocated N pool (gN/m2/s)          
         frootn_to_litter      =>    cnveg_nitrogenflux_inst%frootn_to_litter_patch    , & ! Output: [real(r8) (:) ]  fine root N litterfall (gN/m2/s)                  
         leafc_to_litter_fun   =>    cnveg_carbonflux_inst%leafc_to_litter_fun_patch   , & ! Output:  [real(r8) (:) ]  leaf C litterfall used by FUN (gC/m2/s)
         leafcn_offset         =>    cnveg_state_inst%leafcn_offset_patch               & ! Output:  [real(r8) (:) ]  Leaf C:N used by FUN
         )

      ! The litterfall transfer rate starts at 0.0 and increases linearly
      ! over time, with displayed growth going to 0.0 on the last day of litterfall
      
      do fp = 1,num_soilp
         p = filter_soilp(fp)

         ! only calculate fluxes during offset period
         if (offset_flag(p) == 1._r8) then

            if (offset_counter(p) == dt) then
               t1 = 1.0_r8 / dt
               frootc_to_litter(p) = t1 * frootc(p) + cpool_to_frootc(p)
               
               ! biofuel_harvfrac is only non-zero for prognostic crops.
               leafc_to_litter(p)  = t1 * leafc(p)*(1._r8-biofuel_harvfrac(ivt(p)))  + cpool_to_leafc(p)

               ! this assumes that offset_counter == dt for crops
               ! if this were ever changed, we'd need to add code to the "else"
               if (ivt(p) >= npcropmin) then
                  ! Replenish the seed deficits from grain, if there is enough
                  ! available grain. (If there is not enough available grain, the seed
                  ! deficits will accumulate until there is eventually enough grain to
                  ! replenish them.)
                  grainc_to_seed(p) = t1 * min(-cropseedc_deficit(p), grainc(p))
                  grainn_to_seed(p) = t1 * min(-cropseedn_deficit(p), grainn(p))
                  ! Send the remaining grain to the food product pool
                  grainc_to_food(p) = t1 * grainc(p)  + cpool_to_grainc(p) - grainc_to_seed(p)
                  grainn_to_food(p) = t1 * grainn(p)  + npool_to_grainn(p) - grainn_to_seed(p)
                  
                  ! Cut a certain fraction (i.e., biofuel_harvfrac(ivt(p))) (e.g., biofuel_harvfrac(ivt(p)=70% for bioenergy crops) of leaf C
                  ! and move this fration of leaf C to biofuel C, rather than move it to litter
                  leafc_to_biofuelc(p) = t1 * leafc(p) * biofuel_harvfrac(ivt(p))
                  leafn_to_biofueln(p) = t1 * leafn(p) * biofuel_harvfrac(ivt(p))

                  ! Cut a certain fraction (i.e., biofuel_harvfrac(ivt(p))) (e.g., biofuel_harvfrac(ivt(p)=70% for bioenergy crops) of livestem C
                  ! and move this fration of leaf C to biofuel C, rather than move it to litter
                  livestemc_to_litter(p)   = t1 * livestemc(p)*(1._r8-biofuel_harvfrac(ivt(p)))  + cpool_to_livestemc(p)
                  livestemc_to_biofuelc(p) = t1 * livestemc(p) * biofuel_harvfrac(ivt(p))
                  livestemn_to_biofueln(p) = t1 * livestemn(p) * biofuel_harvfrac(ivt(p))
               end if
            else
               t1 = dt * 2.0_r8 / (offset_counter(p) * offset_counter(p))
               leafc_to_litter(p)  = prev_leafc_to_litter(p)  + t1*(leafc(p)  - prev_leafc_to_litter(p)*offset_counter(p))
               frootc_to_litter(p) = prev_frootc_to_litter(p) + t1*(frootc(p) - prev_frootc_to_litter(p)*offset_counter(p))

            end if
            
            if ( use_fun ) then
               if(leafc_to_litter(p)*dt.gt.leafc(p))then
                   leafc_to_litter(p) = leafc(p)/dt + cpool_to_leafc(p)
               endif
               if(frootc_to_litter(p)*dt.gt.frootc(p))then
                   frootc_to_litter(p) = frootc(p)/dt + cpool_to_frootc(p)
               endif
            end if
            
            
            if ( use_fun ) then
               leafc_to_litter_fun(p)      =  leafc_to_litter(p)
               leafn_to_retransn(p)        =  paid_retransn_to_npool(p) + free_retransn_to_npool(p)
               if (leafn(p).gt.0._r8) then
                  if (leafn(p)-leafn_to_retransn(p)*dt.gt.0._r8) then
                      leafcn_offset(p)     =  leafc(p)/(leafn(p)-leafn_to_retransn(p)*dt)
                  else
                      leafcn_offset(p)     =  leafc(p)/leafn(p)
                  end if
               else
                  leafcn_offset(p)         =  leafcn(ivt(p))
               end if
               leafn_to_litter(p)          =  leafc_to_litter(p)/leafcn_offset(p) - leafn_to_retransn(p)
               leafn_to_litter(p)          =  max(leafn_to_litter(p),0._r8)
               
               denom = ( leafn_to_retransn(p) + leafn_to_litter(p) )
               if ( denom /= 0.0_r8 ) then
                  fr_leafn_to_litter =  leafn_to_litter(p) / ( leafn_to_retransn(p) + leafn_to_litter(p) )
               else if ( leafn_to_litter(p) == 0.0_r8 ) then
                  fr_leafn_to_litter =  0.0_r8
               else
                  fr_leafn_to_litter =  1.0_r8
               end if

            else
               if (CNratio_floating .eqv. .true.) then    
                  fr_leafn_to_litter = 0.5_r8    ! assuming 50% of nitrogen turnover goes to litter
               end if
               ! calculate the leaf N litterfall and retranslocation
               leafn_to_litter(p)   = leafc_to_litter(p)  / lflitcn(ivt(p))
               leafn_to_retransn(p) = (leafc_to_litter(p) / leafcn(ivt(p))) - leafn_to_litter(p)

            end if    

            ! calculate fine root N litterfall (no retranslocation of fine root N)
            frootn_to_litter(p) = frootc_to_litter(p) / frootcn(ivt(p))
            
            if (CNratio_floating .eqv. .true.) then    
               if (leafc(p) == 0.0_r8) then    
                  ntovr_leaf = 0.0_r8    
               else    
                  ntovr_leaf = leafc_to_litter(p) * (leafn(p) / leafc(p))   
               end if   
           
               leafn_to_litter(p)   = fr_leafn_to_litter * ntovr_leaf
               leafn_to_retransn(p) = ntovr_leaf - leafn_to_litter(p)
               if (frootc(p) == 0.0_r8) then    
                   frootn_to_litter(p) = 0.0_r8    
                else    
                   frootn_to_litter(p) = frootc_to_litter(p) * (frootn(p) / frootc(p))   
                end if   
            end if  
            
            if ( use_fun ) then
               if(frootn_to_litter(p)*dt.gt.frootn(p))then
                   frootn_to_litter(p) = frootn(p)/dt
               endif    
            end if

            if (ivt(p) >= npcropmin) then
               ! NOTE(slevis, 2014-12) results in -ve livestemn and -ve totpftn
               !X! livestemn_to_litter(p) = livestemc_to_litter(p) / livewdcn(ivt(p))
               ! NOTE(slevis, 2014-12) Beth Drewniak suggested this instead
               livestemn_to_litter(p) = livestemn(p) / dt * (1 - biofuel_harvfrac(ivt(p)))
            end if

            ! save the current litterfall fluxes
            prev_leafc_to_litter(p)  = leafc_to_litter(p)
            prev_frootc_to_litter(p) = frootc_to_litter(p)
                
         end if ! end if offset period

      end do ! end patch loop

    end associate 

  end subroutine CNOffsetLitterfall

  !-----------------------------------------------------------------------
  subroutine CNBackgroundLitterfall (num_soilp, filter_soilp, &
       cnveg_state_inst, cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! Determines the flux of C and N from displayed pools to litter
    ! pools as the result of background litter fall.
    !
    ! !USES:
    use CNSharedParamsMod   , only : use_fun
    use clm_varctl          , only : CNratio_floating    
    ! !ARGUMENTS:
    implicit none
    integer                       , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                       , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(cnveg_state_type)        , intent(inout) :: cnveg_state_inst
    type(cnveg_carbonstate_type)  , intent(in)    :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type), intent(in)    :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)   , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type) , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    integer :: p            ! indices
    integer :: fp           ! lake filter patch index
    real(r8) :: fr_leafn_to_litter ! fraction of the nitrogen turnover that goes to litter; remaining fraction is retranslocated
    real(r8) :: ntovr_leaf  
    real(r8) :: denom       
    !-----------------------------------------------------------------------

    associate(                                                                     & 
         ivt               =>    patch%itype                                     , & ! Input:  [integer  (:) ]  patch vegetation type                                

         leafcn            =>    pftcon%leafcn                                   , & ! Input:  leaf C:N (gC/gN)                                  
         lflitcn           =>    pftcon%lflitcn                                  , & ! Input:  leaf litter C:N (gC/gN)                           
         frootcn           =>    pftcon%frootcn                                  , & ! Input:  fine root C:N (gC/gN)                             

         bglfr             =>    cnveg_state_inst%bglfr_patch                    , & ! Input:  [real(r8) (:) ]  background litterfall rate (1/s)                  

         leafc             =>    cnveg_carbonstate_inst%leafc_patch              , & ! Input:  [real(r8) (:) ]  (gC/m2) leaf C                                    
         frootc            =>    cnveg_carbonstate_inst%frootc_patch             , & ! Input:  [real(r8) (:) ]  (gC/m2) fine root C                               
         
         leafc_to_litter   =>    cnveg_carbonflux_inst%leafc_to_litter_patch     , & ! Output: [real(r8) (:) ]                                                    
         frootc_to_litter  =>    cnveg_carbonflux_inst%frootc_to_litter_patch    , & ! Output: [real(r8) (:) ]                                                    

         leafn             =>    cnveg_nitrogenstate_inst%leafn_patch            , & ! Input:  [real(r8) (:) ]  (gN/m2) leaf N  
         frootn            =>    cnveg_nitrogenstate_inst%frootn_patch           , & ! Input:  [real(r8) (:) ]  (gN/m2) fine root N 
         leafn_to_litter   =>    cnveg_nitrogenflux_inst%leafn_to_litter_patch   , & ! Output: [real(r8) (:) ]                                                    
         leafn_to_retransn =>    cnveg_nitrogenflux_inst%leafn_to_retransn_patch , & ! Output: [real(r8) (:) ]                                                    
         frootn_to_litter  =>    cnveg_nitrogenflux_inst%frootn_to_litter_patch  , & ! Output: [real(r8) (:) ]                                                    
         leafc_to_litter_fun   => cnveg_carbonflux_inst%leafc_to_litter_fun_patch, & ! Output:  [real(r8) (:) ] leaf C litterfall used by FUN (gC/m2/s)
         leafcn_offset         => cnveg_state_inst%leafcn_offset_patch           , & ! Output:  [real(r8) (:) ] Leaf C:N used by FUN
         free_retransn_to_npool=>    cnveg_nitrogenflux_inst%free_retransn_to_npool_patch  , & ! Input: [real(r8) (:) ] free leaf N to retranslocated N pool (gN/m2/s)          
         paid_retransn_to_npool=>    cnveg_nitrogenflux_inst%retransn_to_npool_patch   & ! Input: [real(r8) (:) ] free leaf N to retranslocated N pool (gN/m2/s)          
         )

      ! patch loop
      do fp = 1,num_soilp
         p = filter_soilp(fp)

         ! only calculate these fluxes if the background litterfall rate is non-zero
         if (bglfr(p) > 0._r8) then
            ! units for bglfr are already 1/s
            leafc_to_litter(p)  = bglfr(p) * leafc(p)
            frootc_to_litter(p) = bglfr(p) * frootc(p)
            if ( use_fun ) then
               leafc_to_litter_fun(p)     = leafc_to_litter(p)
               leafn_to_retransn(p)       = paid_retransn_to_npool(p) + free_retransn_to_npool(p)
               if (leafn(p).gt.0._r8) then
                  if (leafn(p)-leafn_to_retransn(p)*dt.gt.0._r8) then
                     leafcn_offset(p)     = leafc(p)/(leafn(p)-leafn_to_retransn(p)*dt)
                  else
                     leafcn_offset(p)     = leafc(p)/leafn(p)
                  end if
               else
                  leafcn_offset(p)        = leafcn(ivt(p))
               end if
               leafn_to_litter(p)         = leafc_to_litter(p)/leafcn_offset(p) - leafn_to_retransn(p)
               leafn_to_litter(p)         = max(leafn_to_litter(p),0._r8)

               denom = ( leafn_to_retransn(p) + leafn_to_litter(p) )
               if ( denom /= 0.0_r8 ) then
                  fr_leafn_to_litter =  leafn_to_litter(p) / ( leafn_to_retransn(p) + leafn_to_litter(p) )
               else if ( leafn_to_litter(p) == 0.0_r8 ) then
                  fr_leafn_to_litter =  0.0_r8
               else
                  fr_leafn_to_litter =  1.0_r8
               end if


            else
               if (CNratio_floating .eqv. .true.) then    
                  fr_leafn_to_litter = 0.5_r8    ! assuming 50% of nitrogen turnover goes to litter
               end if
               ! calculate the leaf N litterfall and retranslocation
               leafn_to_litter(p)   = leafc_to_litter(p)  / lflitcn(ivt(p))
               leafn_to_retransn(p) = (leafc_to_litter(p) / leafcn(ivt(p))) - leafn_to_litter(p)

            end if    

            ! calculate fine root N litterfall (no retranslocation of fine root N)
            frootn_to_litter(p) = frootc_to_litter(p) / frootcn(ivt(p))
            
            if (CNratio_floating .eqv. .true.) then    
               if (leafc(p) == 0.0_r8) then    
                  ntovr_leaf = 0.0_r8    
               else    
                  ntovr_leaf = leafc_to_litter(p) * (leafn(p) / leafc(p))   
               end if   
           
               leafn_to_litter(p)   = fr_leafn_to_litter * ntovr_leaf
               leafn_to_retransn(p) = ntovr_leaf - leafn_to_litter(p)
               if (frootc(p) == 0.0_r8) then    
                   frootn_to_litter(p) = 0.0_r8    
                else    
                   frootn_to_litter(p) = frootc_to_litter(p) * (frootn(p) / frootc(p))   
                end if   
            end if    

            if ( use_fun ) then
               if(frootn_to_litter(p)*dt.gt.frootn(p))then
                    frootn_to_litter(p) = frootn(p)/dt
               endif
            end if

         end if

      end do

    end associate 

  end subroutine CNBackgroundLitterfall

  !-----------------------------------------------------------------------
  subroutine CNLivewoodTurnover (num_soilp, filter_soilp, &
       cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! Determines the flux of C and N from live wood to
    ! dead wood pools, for stem and coarse root.
    !
    use CNSharedParamsMod, only: use_fun
    use clm_varctl          , only : CNratio_floating    
    ! !ARGUMENTS:
    integer                        , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                        , intent(in)    :: filter_soilp(:) ! filter for soil patches
    type(cnveg_carbonstate_type)   , intent(in)    :: cnveg_carbonstate_inst
    type(cnveg_nitrogenstate_type) , intent(in)    :: cnveg_nitrogenstate_inst
    type(cnveg_carbonflux_type)    , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type)  , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    integer :: p            ! indices
    integer :: fp           ! lake filter patch index
    real(r8):: ctovr        ! temporary variable for carbon turnover
    real(r8):: ntovr        ! temporary variable for nitrogen turnover
    !-----------------------------------------------------------------------

    associate(                                                                                   & 
         ivt                      =>    patch%itype                                              , & ! Input:  [integer  (:) ]  patch vegetation type                                

         woody                    =>    pftcon%woody                                           , & ! Input:  binary flag for woody lifeform (1=woody, 0=not woody)
         livewdcn                 =>    pftcon%livewdcn                                        , & ! Input:  live wood (phloem and ray parenchyma) C:N (gC/gN) 
         deadwdcn                 =>    pftcon%deadwdcn                                        , & ! Input:  dead wood (xylem and heartwood) C:N (gC/gN)       

         livestemc                =>    cnveg_carbonstate_inst%livestemc_patch                 , & ! Input:  [real(r8) (:) ]  (gC/m2) live stem C                               
         livecrootc               =>    cnveg_carbonstate_inst%livecrootc_patch                , & ! Input:  [real(r8) (:) ]  (gC/m2) live coarse root C                        

         livestemn                =>    cnveg_nitrogenstate_inst%livestemn_patch               , & ! Input:  [real(r8) (:) ]  (gN/m2) live stem N                               
         livecrootn               =>    cnveg_nitrogenstate_inst%livecrootn_patch              , & ! Input:  [real(r8) (:) ]  (gN/m2) live coarse root N                        
         
         livestemc_to_deadstemc   =>    cnveg_carbonflux_inst%livestemc_to_deadstemc_patch     , & ! Output: [real(r8) (:) ]                                                    
         livecrootc_to_deadcrootc =>    cnveg_carbonflux_inst%livecrootc_to_deadcrootc_patch   , & ! Output: [real(r8) (:) ]                                                    

         livestemn_to_deadstemn   =>    cnveg_nitrogenflux_inst%livestemn_to_deadstemn_patch   , & ! Output: [real(r8) (:) ]                                                    
         livestemn_to_retransn    =>    cnveg_nitrogenflux_inst%livestemn_to_retransn_patch    , & ! Output: [real(r8) (:) ]                                                    
         livecrootn_to_deadcrootn =>    cnveg_nitrogenflux_inst%livecrootn_to_deadcrootn_patch , & ! Output: [real(r8) (:) ]                                                    
         livecrootn_to_retransn   =>    cnveg_nitrogenflux_inst%livecrootn_to_retransn_patch     & ! Output: [real(r8) (:) ]                                                    
         )



      ! patch loop
      do fp = 1,num_soilp
         p = filter_soilp(fp)

         ! only calculate these fluxes for woody types
         if (woody(ivt(p)) > 0._r8) then

            ! live stem to dead stem turnover

            ctovr = livestemc(p) * lwtop
            ntovr = ctovr / livewdcn(ivt(p))
            livestemc_to_deadstemc(p) = ctovr
            livestemn_to_deadstemn(p) = ctovr / deadwdcn(ivt(p))
            
            if (CNratio_floating .eqv. .true.) then    
               if (livestemc(p) == 0.0_r8) then    
                   ntovr = 0.0_r8    
                else    
                   ntovr = ctovr * (livestemn(p) / livestemc(p))   
                end if   

                livestemn_to_deadstemn(p) = 0.5_r8 * ntovr   ! assuming 50% goes to deadstemn 
            end if    
            
            livestemn_to_retransn(p)  = ntovr - livestemn_to_deadstemn(p)

            ! live coarse root to dead coarse root turnover

            ctovr = livecrootc(p) * lwtop
            ntovr = ctovr / livewdcn(ivt(p))
            livecrootc_to_deadcrootc(p) = ctovr
            livecrootn_to_deadcrootn(p) = ctovr / deadwdcn(ivt(p))
            
            if (CNratio_floating .eqv. .true.) then    
              if (livecrootc(p) == 0.0_r8) then    
                  ntovr = 0.0_r8    
               else    
                  ntovr = ctovr * (livecrootn(p) / livecrootc(p))   
               end if   

               livecrootn_to_deadcrootn(p) = 0.5_r8 * ntovr   ! assuming 50% goes to deadstemn 
            end if    
            
            livecrootn_to_retransn(p)  = ntovr - livecrootn_to_deadcrootn(p)
            if(use_fun)then
               !TURNED OFF FLUXES TO CORRECT N ACCUMULATION ISSUE. RF. Oct 2015. 
               livecrootn_to_retransn(p) = 0.0_r8
               livestemn_to_retransn(p)  = 0.0_r8
            endif

         end if

      end do

    end associate 

  end subroutine CNLivewoodTurnover

  !-----------------------------------------------------------------------
  subroutine CNCropHarvestToProductPools(bounds, num_soilp, filter_soilp, num_soilc, filter_soilc, &
       cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
    !
    ! !DESCRIPTION:
    ! If using prognostic crop, then move any necessary harvested amounts into fluxes
    ! destined for the product pools.
    !
    ! !USES:
    use clm_varctl    , only : use_crop
    use clm_varctl    , only : use_grainproduct
    use subgridAveMod , only : p2c
    !
    ! !ARGUMENTS:
    type(bounds_type)             , intent(in)    :: bounds
    integer                       , intent(in)    :: num_soilp       ! number of soil patches in filter
    integer                       , intent(in)    :: filter_soilp(:) ! filter for soil patches
    integer                       , intent(in)    :: num_soilc       ! number of soil columns in filter
    integer                       , intent(in)    :: filter_soilc(:) ! filter for soil columns
    type(cnveg_carbonflux_type)   , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type) , intent(inout) :: cnveg_nitrogenflux_inst
    !
    ! !LOCAL VARIABLES:
    integer :: fp, p

    character(len=*), parameter :: subname = 'CNCropHarvestToProductPools'
    !-----------------------------------------------------------------------

    if (use_crop) then
       do fp = 1, num_soilp
          p = filter_soilp(fp)
          cnveg_carbonflux_inst%grainc_to_cropprodc_patch(p) = cnveg_carbonflux_inst%leafc_to_biofuelc_patch(p) + &
               cnveg_carbonflux_inst%livestemc_to_biofuelc_patch(p)
          cnveg_nitrogenflux_inst%grainn_to_cropprodn_patch(p) = cnveg_nitrogenflux_inst%leafn_to_biofueln_patch(p) + &
               cnveg_nitrogenflux_inst%livestemn_to_biofueln_patch(p)
       end do

       if (use_grainproduct) then
          do fp = 1, num_soilp
             p = filter_soilp(fp)
             cnveg_carbonflux_inst%grainc_to_cropprodc_patch(p) = cnveg_carbonflux_inst%grainc_to_cropprodc_patch(p) + &
                  cnveg_carbonflux_inst%grainc_to_food_patch(p)
             cnveg_nitrogenflux_inst%grainn_to_cropprodn_patch(p) = cnveg_nitrogenflux_inst%grainn_to_cropprodn_patch(p) + &
                  cnveg_nitrogenflux_inst%grainn_to_food_patch(p)
          end do
       end if
       
       call p2c (bounds, num_soilc, filter_soilc, &
            cnveg_carbonflux_inst%grainc_to_cropprodc_patch(bounds%begp:bounds%endp), &
            cnveg_carbonflux_inst%grainc_to_cropprodc_col(bounds%begc:bounds%endc))

       call p2c (bounds, num_soilc, filter_soilc, &
            cnveg_nitrogenflux_inst%grainn_to_cropprodn_patch(bounds%begp:bounds%endp), &
            cnveg_nitrogenflux_inst%grainn_to_cropprodn_col(bounds%begc:bounds%endc))

    end if

  end subroutine CNCropHarvestToProductPools

  !-----------------------------------------------------------------------
  subroutine CNLitterToColumn (bounds, num_soilc, filter_soilc,         &
       cnveg_state_inst,cnveg_carbonflux_inst, cnveg_nitrogenflux_inst, &
       leaf_prof_patch, froot_prof_patch)
    !
    ! !DESCRIPTION:
    ! called at the end of cn_phenology to gather all patch-level litterfall fluxes
    ! to the column level and assign them to the three litter pools
    !
    ! !USES:
    use clm_varpar , only : max_patch_per_col, nlevdecomp
    use pftconMod  , only : npcropmin
    use clm_varctl , only : use_grainproduct
    !
    ! !ARGUMENTS:
    type(bounds_type)               , intent(in)    :: bounds
    integer                         , intent(in)    :: num_soilc       ! number of soil columns in filter
    integer                         , intent(in)    :: filter_soilc(:) ! filter for soil columns
    type(cnveg_state_type)          , intent(in)    :: cnveg_state_inst
    type(cnveg_carbonflux_type)     , intent(inout) :: cnveg_carbonflux_inst
    type(cnveg_nitrogenflux_type)   , intent(inout) :: cnveg_nitrogenflux_inst
    real(r8)                        , intent(in)    :: leaf_prof_patch(bounds%begp:,1:)
    real(r8)                        , intent(in)    :: froot_prof_patch(bounds%begp:,1:)
    !
    ! !LOCAL VARIABLES:
    integer :: fc,c,pi,p,j       ! indices
    !-----------------------------------------------------------------------

    SHR_ASSERT_ALL_FL((ubound(leaf_prof_patch)   == (/bounds%endp,nlevdecomp_full/)), sourcefile, __LINE__)
    SHR_ASSERT_ALL_FL((ubound(froot_prof_patch)  == (/bounds%endp,nlevdecomp_full/)), sourcefile, __LINE__)

    associate(                                                                                & 
         leaf_prof                 => leaf_prof_patch                                       , & ! Input:  [real(r8) (:,:) ]  (1/m) profile of leaves                         
         froot_prof                => froot_prof_patch                                      , & ! Input:  [real(r8) (:,:) ]  (1/m) profile of fine roots                     

         ivt                       => patch%itype                                             , & ! Input:  [integer  (:)   ]  patch vegetation type                                
         wtcol                     => patch%wtcol                                             , & ! Input:  [real(r8) (:)   ]  weight (relative to column) for this patch (0-1)    

         lf_flab                   => pftcon%lf_flab                                        , & ! Input:  leaf litter labile fraction                       
         lf_fcel                   => pftcon%lf_fcel                                        , & ! Input:  leaf litter cellulose fraction                    
         lf_flig                   => pftcon%lf_flig                                        , & ! Input:  leaf litter lignin fraction                       
         fr_flab                   => pftcon%fr_flab                                        , & ! Input:  fine root litter labile fraction                  
         fr_fcel                   => pftcon%fr_fcel                                        , & ! Input:  fine root litter cellulose fraction               
         fr_flig                   => pftcon%fr_flig                                        , & ! Input:  fine root litter lignin fraction                  

         leafc_to_litter           => cnveg_carbonflux_inst%leafc_to_litter_patch           , & ! Input:  [real(r8) (:)   ]  leaf C litterfall (gC/m2/s)                       
         frootc_to_litter          => cnveg_carbonflux_inst%frootc_to_litter_patch          , & ! Input:  [real(r8) (:)   ]  fine root N litterfall (gN/m2/s)                  
         livestemc_to_litter       => cnveg_carbonflux_inst%livestemc_to_litter_patch       , & ! Input:  [real(r8) (:)   ]  live stem C litterfall (gC/m2/s)                  
         grainc_to_food            => cnveg_carbonflux_inst%grainc_to_food_patch            , & ! Input:  [real(r8) (:)   ]  grain C to food (gC/m2/s)                            
         phenology_c_to_litr_met_c => cnveg_carbonflux_inst%phenology_c_to_litr_met_c_col   , & ! Output: [real(r8) (:,:) ]  C fluxes associated with phenology (litterfall and crop) to litter metabolic pool (gC/m3/s)
         phenology_c_to_litr_cel_c => cnveg_carbonflux_inst%phenology_c_to_litr_cel_c_col   , & ! Output: [real(r8) (:,:) ]  C fluxes associated with phenology (litterfall and crop) to litter cellulose pool (gC/m3/s)
         phenology_c_to_litr_lig_c => cnveg_carbonflux_inst%phenology_c_to_litr_lig_c_col   , & ! Output: [real(r8) (:,:) ]  C fluxes associated with phenology (litterfall and crop) to litter lignin pool (gC/m3/s)

         livestemn_to_litter       => cnveg_nitrogenflux_inst%livestemn_to_litter_patch     , & ! Input:  [real(r8) (:)   ]  livestem N to litter (gN/m2/s)                    
         grainn_to_food            => cnveg_nitrogenflux_inst%grainn_to_food_patch          , & ! Input:  [real(r8) (:)   ]  grain N to food (gN/m2/s) 
         leafn_to_litter           => cnveg_nitrogenflux_inst%leafn_to_litter_patch         , & ! Input:  [real(r8) (:)   ]  leaf N litterfall (gN/m2/s)                       
         frootn_to_litter          => cnveg_nitrogenflux_inst%frootn_to_litter_patch        , & ! Input:  [real(r8) (:)   ]  fine root N litterfall (gN/m2/s)                  
         phenology_n_to_litr_met_n => cnveg_nitrogenflux_inst%phenology_n_to_litr_met_n_col , & ! Output: [real(r8) (:,:) ]  N fluxes associated with phenology (litterfall and crop) to litter metabolic pool (gN/m3/s)
         phenology_n_to_litr_cel_n => cnveg_nitrogenflux_inst%phenology_n_to_litr_cel_n_col , & ! Output: [real(r8) (:,:) ]  N fluxes associated with phenology (litterfall and crop) to litter cellulose pool (gN/m3/s)
         phenology_n_to_litr_lig_n => cnveg_nitrogenflux_inst%phenology_n_to_litr_lig_n_col   & ! Output: [real(r8) (:,:) ]  N fluxes associated with phenology (litterfall and crop) to litter lignin pool (gN/m3/s)
         )
    
      do j = 1, nlevdecomp
         do pi = 1,max_patch_per_col
            do fc = 1,num_soilc
               c = filter_soilc(fc)

               if ( pi <=  col%npatches(c) ) then
                  p = col%patchi(c) + pi - 1
                  if (patch%active(p)) then

                     ! leaf litter carbon fluxes
                     phenology_c_to_litr_met_c(c,j) = phenology_c_to_litr_met_c(c,j) &
                          + leafc_to_litter(p) * lf_flab(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                     phenology_c_to_litr_cel_c(c,j) = phenology_c_to_litr_cel_c(c,j) &
                          + leafc_to_litter(p) * lf_fcel(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                     phenology_c_to_litr_lig_c(c,j) = phenology_c_to_litr_lig_c(c,j) &
                          + leafc_to_litter(p) * lf_flig(ivt(p)) * wtcol(p) * leaf_prof(p,j)

                     ! leaf litter nitrogen fluxes
                     phenology_n_to_litr_met_n(c,j) = phenology_n_to_litr_met_n(c,j) &
                          + leafn_to_litter(p) * lf_flab(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                     phenology_n_to_litr_cel_n(c,j) = phenology_n_to_litr_cel_n(c,j) &
                          + leafn_to_litter(p) * lf_fcel(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                     phenology_n_to_litr_lig_n(c,j) = phenology_n_to_litr_lig_n(c,j) &
                          + leafn_to_litter(p) * lf_flig(ivt(p)) * wtcol(p) * leaf_prof(p,j)

                     ! fine root litter carbon fluxes
                     phenology_c_to_litr_met_c(c,j) = phenology_c_to_litr_met_c(c,j) &
                          + frootc_to_litter(p) * fr_flab(ivt(p)) * wtcol(p) * froot_prof(p,j)
                     phenology_c_to_litr_cel_c(c,j) = phenology_c_to_litr_cel_c(c,j) &
                          + frootc_to_litter(p) * fr_fcel(ivt(p)) * wtcol(p) * froot_prof(p,j)
                     phenology_c_to_litr_lig_c(c,j) = phenology_c_to_litr_lig_c(c,j) &
                          + frootc_to_litter(p) * fr_flig(ivt(p)) * wtcol(p) * froot_prof(p,j)

                     ! fine root litter nitrogen fluxes
                     phenology_n_to_litr_met_n(c,j) = phenology_n_to_litr_met_n(c,j) &
                          + frootn_to_litter(p) * fr_flab(ivt(p)) * wtcol(p) * froot_prof(p,j)
                     phenology_n_to_litr_cel_n(c,j) = phenology_n_to_litr_cel_n(c,j) &
                          + frootn_to_litter(p) * fr_fcel(ivt(p)) * wtcol(p) * froot_prof(p,j)
                     phenology_n_to_litr_lig_n(c,j) = phenology_n_to_litr_lig_n(c,j) &
                          + frootn_to_litter(p) * fr_flig(ivt(p)) * wtcol(p) * froot_prof(p,j)

                     ! agroibis puts crop stem litter together with leaf litter
                     ! so I've used the leaf lf_f* parameters instead of making
                     ! new ones for now (slevis)
                     ! also for simplicity I've put "food" into the litter pools

                     if (ivt(p) >= npcropmin) then ! add livestemc to litter
                        ! stem litter carbon fluxes
                        phenology_c_to_litr_met_c(c,j) = phenology_c_to_litr_met_c(c,j) &
                             + livestemc_to_litter(p) * lf_flab(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                        phenology_c_to_litr_cel_c(c,j) = phenology_c_to_litr_cel_c(c,j) &
                             + livestemc_to_litter(p) * lf_fcel(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                        phenology_c_to_litr_lig_c(c,j) = phenology_c_to_litr_lig_c(c,j) &
                             + livestemc_to_litter(p) * lf_flig(ivt(p)) * wtcol(p) * leaf_prof(p,j)

                        ! stem litter nitrogen fluxes
                        phenology_n_to_litr_met_n(c,j) = phenology_n_to_litr_met_n(c,j) &
                             + livestemn_to_litter(p) * lf_flab(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                        phenology_n_to_litr_cel_n(c,j) = phenology_n_to_litr_cel_n(c,j) &
                             + livestemn_to_litter(p) * lf_fcel(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                        phenology_n_to_litr_lig_n(c,j) = phenology_n_to_litr_lig_n(c,j) &
                             + livestemn_to_litter(p) * lf_flig(ivt(p)) * wtcol(p) * leaf_prof(p,j)

                        if (.not. use_grainproduct) then
                         ! grain litter carbon fluxes
                         phenology_c_to_litr_met_c(c,j) = phenology_c_to_litr_met_c(c,j) &
                              + grainc_to_food(p) * lf_flab(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                         phenology_c_to_litr_cel_c(c,j) = phenology_c_to_litr_cel_c(c,j) &
                              + grainc_to_food(p) * lf_fcel(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                         phenology_c_to_litr_lig_c(c,j) = phenology_c_to_litr_lig_c(c,j) &
                              + grainc_to_food(p) * lf_flig(ivt(p)) * wtcol(p) * leaf_prof(p,j)
 
                         ! grain litter nitrogen fluxes
                         phenology_n_to_litr_met_n(c,j) = phenology_n_to_litr_met_n(c,j) &
                              + grainn_to_food(p) * lf_flab(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                         phenology_n_to_litr_cel_n(c,j) = phenology_n_to_litr_cel_n(c,j) &
                              + grainn_to_food(p) * lf_fcel(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                         phenology_n_to_litr_lig_n(c,j) = phenology_n_to_litr_lig_n(c,j) &
                              + grainn_to_food(p) * lf_flig(ivt(p)) * wtcol(p) * leaf_prof(p,j)
                        end if


                     end if
                  end if
               end if

            end do

         end do
      end do

    end associate 

  end subroutine CNLitterToColumn

end module CNPhenologyMod
